{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Coarsen in Two Dimensions\n",
    "\n",
    "We describe the coarsening algorithm implemented in `coarsen.m` for two dimensional bisection grids. We assume when bisecting a triangle, the left child is stored in a prori of the right child. More details can be found in the paper,\n",
    "\n",
    "- L. Chen and C-S. Zhang. [A coarsening algorithm on adaptive grids by newest vertex bisection and its applications](http://math.uci.edu/~chenlong/CZ2007.html). *Journal of Computational Mathematics*, 28(6):767-789, 2010."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compatible Bisection\n",
    "Every bisection grids can be obtained from an initial grid with a\n",
    "sequence of compatible bisections. In two dimensions, there are two types\n",
    "of compatible bisections.\n",
    "\n",
    "- Case 1: A pair of triangles are divided into four triangles by adding the middle point their common edge."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBgcVeqVoUgAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNDowNzoyMVeWFJ0AABAN\nSURBVHic7do/q6XXecbhZ0I6g1Eb8MCeLuBA5MZmVJhRp26iRukk9AlikL6Aepv4G4i4D7Y6g4tR\nZWFV+QgbJkXKYcClmRTHHB+fP/vsP+/7rnutdV3V4VSreG5+e8+cJ+/evSsAaO0fWj8AAKoECYAQ\nggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIA\nEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAk\nACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIgg\nSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQ\nQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEEC\nIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiC\nBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgAR\nBAmACIIEQARBAoay3+9bP4EzCdJc9vv9559/brGMar/fP3v27NmzZ468R4I0kf1+/5Of/+ub9//v\n409fmivjubrw3fOn7/30B468R4I0i6ut/tt/fvT+Jz9+8cufmSuDub7w3Qc/eu/pDx15jwRpCje2\n+rSqzJXB3LrwcuR9EqTx3d1qmSsDuffCy5F3SJAG99BWy1wZwoELL0feG0Ea2eGtlrnSuUcvvBx5\nVwRpWMdstcyVbh154eXI+yFIYzp+q2WudOikCy9H3glBGtCpWy1zpStnXHg58h4I0mjO22qZK504\n+8LLkccTpKFcstUyV+JdeOHlyLMJ0jgu32qZK8EWufBy5MEEaRBLbbXMlUgLXng58lSCNIJlt1rm\nSpjFL7wceSRB6t4aWy1zJcZKF16OPI8g9W29rZa5EmDVCy9HHkaQOrb2VstcaWqDCy9HnkSQerXN\nVstcaWSzCy9HHkOQurTlVstc2dzGF16OPIMg9Wf7rZa5sqEmF16OPIAgdabVVstc2UTDCy9H3pog\n9aTtVstcWVnzCy9H3pQgdSNhq2WurCbkwsuRtyNIfcjZapkrK4i68HLkjQhSB9K2WubKogIvvBx5\nC4KULnOrZa4sJPbCy5FvTpCiJW+1zJWLhV94OfJtCVKu/K2WuXKBLi68HPmGBClUL1stc+UsHV14\nOfKtCFKivrZa5sqJurvwcuSbEKQ4PW61zJWjdXrh5cjXJ0hZ+t1qmStH6PrCy5GvTJCC9L7VMlcO\nGuDCy5GvSZBSjLHVMlceMMyFlyNfjSBFGGmrZa7cMdiFlyNfhyC1N95Wy1y5YcgLL0e+AkFqbNSt\nlrlSVUNfeDnypQlSS2Nvtcx1esNfeDnyRQlSMzNstcx1YpNceDny5QhSG/Nstcx1SlNdeDnyhQhS\nA7Nttcx1MhNeeDnyJQjS1ubcapnrNKa98HLkFxOkTc281TLXCUx+4eXILyNI27HVMtehufArjvxs\ngrQRW71mrkNy4Tc58vMI0hZs9RZzHYwLv8uRn0GQVmer9zLXYbjwhzjyUwnSumz1AHMdgAs/zJGf\nRJBWZKuPMteuufBjOPLjCdJabPVI5topF348R34kQVqFrZ7EXLvjwk/lyI8hSMuz1TOYa0dc+Hkc\n+aMEaWG2ejZz7YILv4QjP0yQlmSrFzLXcC78co78AEFajK0uwlxjufClOPKHCNIybHVB5hrIhS/L\nkd9LkBZgq4sz1ygufA2O/C5BupStrsRcQ7jw9TjyWwTpIra6KnNtzoWvzZHfJEjns9UNmGtDLnwb\njvyaIJ3JVjdjrk248C058iuCdA5b3Zi5bsyFb8+RlyCdwVabMNfNuPBWHLkgncZWGzLXDbjwtiY/\nckE6ga02N/lc1+bCE8x85IJ0LFsNMfNcV+XCc0x75IJ0FFuNMu1c1+PC08x55IL0OFsNNOdcV+LC\nM0145IL0CFuNNeFc1+DCk8125IJ0iK2Gm22ui3Ph+aY6ckF6kK12Yaq5LsuF92KeIxek+9lqR+aZ\n64JceF8mOXJBuoetdmeSuS7FhfdohiMXpNtstVMzzHURLrxfwx+5IP0dW+3a8HO9nAvv3dhHLkh/\nY6sDGHuuF3LhYxj4yAXpr2x1GAPP9RIufCSjHrkgVdnqcEad69lc+HiGPHJBstUxDTnX87jwUY13\n5LMHyVYHNt5cz+DCxzbYkU8dJFsd3mBzPZULn8FIRz5vkGx1EiPN9SQufB7DHPmkQbLVqQwz1+O5\n8NmMceQzBslWJzTGXI/kwuc0wJFPFyRbndYAcz2GC59Z70c+V5BsdXK9z/VRLpyuj3yiINkq1flc\nD3PhXOn3yGcJkq1yrd+5HuDCuanTI58iSLbKLZ3O9SEunLt6PPLxg2Sr3KvHud7Lhd+ya/2AHN0d\n+ZN37961fsOKbJXD3rx+++2Xf/rtb77Z7Xat33IOF37Te6/fvvzF73ffva6q/fOn3/z6ozdPf9j6\nUe11dOQjf0OyVR7V3UfIm1z4Ldc1qqrdd69f/uL3bd8ToqMjH/kb0pMnT3bPn+4++FHrh5Duzeu3\nb77/86tXr/I/Qt7kwm959avvbv3mwy+eN3lJoC6OfORvSLvd7j1f2BmXC79pd/QvifWPrR+wolev\nXn386cv3v/gXo+WAN6/f/u7f//D1118nf3K8lwu/Zf/H/73+J7uq2j9/uvvig1279+To5cgH/4b0\n29988+2Xf3rz+m3rtxDqeqgvXrxo/ZaTufBbvvn1R//zyY+vft4/f/pf//1J2/eE6OjIR/4/pCv7\n/f7jT1+++OXPfIrklo6GeoALv2X/qz/uq1588UHrh0To68hH/oZ0xadI7tXXUA9w4bfsWz8gR3dH\nPn6QymK5o7uhHubCuavHI58iSGWx3NDjUB/lwrmp0yOfJUhlsVRVt0M9hgvnSr9HPlGQymKn1+9Q\nj+TC6frI5wpSWezEuh7q8Vz4zHo/8umCVBY7pd6HehIXPqcBjnzGIJXFTmaAoZ7Khc9mjCOfNEhl\nsdMYY6hncOHzGObI5w1SWewEhhnqeVz4DEY68qmDVBY7tJGGejYXPrbBjnz2IJXFDmqwoV7ChY9q\nvCMXpCqLHc54Q72QCx/PkEcuSH9lscMYcqiXc+EjGfXIBelvLHYAow51ES58DAMfuSD9HYvt2sBD\nXYoL793YRy5It1lsp8Ye6oJceL+GP3JBuofFdmf4oS7LhfdohiMXpPtZbEdmGOriXHhfJjlyQXqQ\nxXZhkqGuwYX3Yp4jF6RDLDbcPENdiQvPN9WRC9IjLDbWVENdjwtPNtuRC9LjLDbQbENdlQvPNOGR\nC9JRLDbKhENdmwtPM+eRC9KxLDbEnEPdgAvPMe2RC9IJLLa5aYe6DReeYOYjF6TTWGxDMw91My68\nrcmPXJBOZrFNTD7ULbnwVhy5IJ3DYjdmqBtz4dtz5CVIZ7PYzRhqEy58S478iiCdz2I3YKgNufBt\nOPJrgnQRi12VoTbnwtfmyG8SpEtZ7EoMNYQLX48jv0WQFmCxizPUKC58DY78LkFahsUuyFADufBl\nOfJ7CdJiLHYRhhrLhS/FkT9EkJZksRcy1HAu/HKO/ABBWpjFns1Qu+DCL+HIDxOk5VnsGQy1Iy78\nPI78UYK0Cos9iaF2x4WfypEfQ5DWYrFHMtROufDjOfIjCdKKLPZRhto1F34MR348QVqXxR5gqANw\n4Yc58pMI0uos9l6GOgwX/hBHfipB2oLF3mKog3HhdznyMwjSRiz2mqEOyYXf5MjPI0jbsdgy1KG5\n8CuO/GyCtKnJF2uow5v8wsuRX0aQtjbtYg11EtNeeDnyiwlSAxMu1lCnMuGFlyNfgiC1MdViDXVC\nU114OfKFCFIzkyzWUKc1yYWXI1+OILU0/GINdXLDX3g58kUJUmMDL9ZQqaEvvBz50gSpvSEXa6hc\nG/LCy5GvQJAiDLZYQ+WWwS68HPk6BCnFMIs1VO41zIWXI1+NIAUZYLGGygEDXHg58jUJUpauF2uo\nPKrrCy9HvjJBitPpYg2VI3V64eXI1ydIibpbrKFyku4uvBz5JgQpVEeLNVTO0NGFlyPfiiDl6mKx\nhsrZurjwcuQbEqRo4Ys1VC4UfuHlyLclSOliF2uoLCL2wsuRb06QOhC4WENlQYEXXo68BUHqQ9Ri\nDZXFRV14OfJGBKkbIYs1VFYScuHlyNsRpJ40X6yhsqrmF16OvClB6kzDxRoqG2jbJEfeliD1p8li\nDZXNtGqSI29OkLq08WINlY1t3yRHnkCQerXZYg2VJrZskiMPIUgd22CxhkpD2zTJkecQpL6tulhD\npbm1m+TIowhS91ZarKESYr0mOfI0gjSCxRdrqERZo0mOPJAgDWLBxRoqgZZtkiPPJEjjWGSxhkqs\npZrkyGMJ0lAuXKyhEu7yJjnyZII0mrMXa6h04ZImOfJwgjSgMxZrqHTkvCY58nyCNKaTFmuodOfU\nJjnyLgjSsI5crKHSqeOb5Mh7IUgje3SxhkrXjmmSI++IIA3uwGINlQEcbpIj74sgje/exRoqw3io\nSY68O4I0hVuLNVQGc7dJjrxHT969e9f6DWxkv99//OnL9//jn7/98ntDZTxXF17/9JeqevP9nx15\ndwRpLvv9/quvvvrss88MlSHt9/sPP/ywqtSoR4IEDGW/31fVbrdr/A5OJ0gARPBHDQBEECQAIggS\nABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQ\nJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACI\nIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiAB\nEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARB\nAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAI\nggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIA\nEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAk\nACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIgg\nSABE+H9E4Z5dq2WzpAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "node = [1,0; 0,1; -1,0; 0,-1];\n",
    "elem = [2 3 1; 4 1 3];\n",
    "figure(1); \n",
    "subplot(1,2,1); showmesh(node,elem);\n",
    "[node,elem] = bisect(node,elem,[1 2]);\n",
    "subplot(1,2,2); showmesh(node,elem); \n",
    "findnode(node,5,'noindex');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Case 2: A triangle containing a boundary edge is divided into two triangles by adding the middle point the boundary edge.\n",
    "\n",
    "The nodes added by compatible bisections are called good-to-carsen nodes\n",
    "or simpley *good* nodes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBgoAotby9AAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNDoxMDowMLIQV4AAAAvF\nSURBVHic7dm9jhxZAYbhGiBDWCPiHdEdgTYBEpA3QLYQYjMvCeFi7sASRCS2cxB7ByuRESC02SZo\nJrMgIoSs0XABgyXiITir2aE9Pz3d9fOdc54n6uroC6r09qk+ury8HABgaV9ZegAADIMgARBCkACI\nIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiAB\nEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARB\nAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAI\nggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIA\nEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAk\nACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIgg\nSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQ\nQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEEC\nIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiC\nBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgAR\nBAmACIIEQARBAiCCIAEQQZCApmw2m6UnsCdB6svZ2dl6vX79+vXSQ2ASm81mvV6v12tZqpEgdeTs\n7OyXv/n4oz/+5Oyff9Ek2rPZbL7/o++uHp8c/+DrP/v4mSZVR5B6UWr0/E8/Pz559L0X39EkGlNq\n9NHvP1x98N7xyaMnv/2hJlVHkLpwVaNyqUk05lqNTso3mlQjQWrfVo0KTaIZ79ao0KTqCFLjbqxR\noUk04LYaFZpUF0Fq2R01KjSJqt1do0KTKiJIzbq3RoUmUaldalRoUi0EqU071qjQJKqze40KTaqC\nIDXoQTUqNImKPLRGhSblE6TW7FGjQpOown41KjQpnCA1Ze8aFZpEuENqVGhSMkFqx4E1KjSJWIfX\nqNCkWILUiFFqVGgSgcaqUaFJmQSpBSPWqNAkooxbo0KTAglS9UavUaFJhJiiRoUmpRGkuk1Uo0KT\nWNx0NSo0KYogVWzSGhWaxIKmrlGhSTkEqVYz1KjQJBYxT40KTQohSFWarUaFJjGzOWtUaFICQarP\nzDUqNInZzF+jQpMWJ0iVWaRGhSYxg6VqVGjSsgSpJgvWqNAkJrVsjQpNWpAgVWPxGhWaxEQSalRo\n0lIEqQ4hNSo0idHl1KjQpEUIUgWialRoEiNKq1GhSfMTpHSBNSo0iVFk1qjQpJkJUrTYGhWaxIGS\na1Ro0pwEKVd4jQpNYm/5NSo0aTaCFKqKGhWaxB5qqVGhSfMQpEQV1ajQJB6krhoVmjQDQYpTXY0K\nTWJHNdao0KSpCVKWSmtUaBL3qrdGhSZNSpCCVF2jQpO4Q+01KjRpOoKUooEaFZrEjdqoUaFJExGk\nCM3UqNAktrRUo0KTpiBIy2usRoUmcaW9GhWaNDpBWliTNSo0iaHdGhWaNC5BWlLDNSo0qXNt16jQ\npBEJ0mKar1GhSd3qoUaFJo1FkJbRSY0KTepQPzUqNGkUgrSArmpUaFJXeqtRoUmHE6S5dVijQpM6\n0WeNCk06kCDNqtsaFZrUvJ5rVGjSIQRpPp3XqNCkhqlRoUl7E6SZqNEVTWqSGl2nSfsRpDmo0RZN\naowavUuT9iBIk1OjG2lSM9ToNpr0UII0LTW6gyY1QI3upkkPIkgTUqN7aVLV1GgXmrQ7QZqKGu1I\nkyqlRrvTpB0J0iTU6EE0qTpq9FCatAtBGp8a7UGTKqJG+9GkewnSyNRob5pUBTU6hCbdTZDGpEYH\n0qRwanQ4TbqDII1GjUahSbHUaCyadBtBGocajUiTAqnRuDTpRoI0AjUanSZFUaMpaNK7BOlQajQR\nTQqhRtPRpC2CdBA1mpQmLU6NpqZJ1wnS/tRoBpq0IDWahyZdEaQ9qdFsNGkRajQnTSoEaR9qNDNN\nmpkazU+TBkHagxotQpNmo0ZL0SRBehg1WpAmzUCNltV5kwTpAdRocZo0KTVK0HOTBGlXahRCkyai\nRjm6bZIg7USNomjS6NQoTZ9NEqT7qVEgTRqRGmXqsEmCdA81iqVJo1CjZL01SZDuokbhNOlAapSv\nqyYJ0q3UqAqatDc1qkU/TRKkm6lRRTRpD2pUl06aJEg3UKPqaNKDqFGNemiSIG1To0pp0o7UqF7N\nN0mQ/o8aVU2T7qVGtWu7SYL0JTVqgCbdQY3a0HCTBOkLatQMTbqRGrWk1SYJ0jCoUXM0aYsatafJ\nJgmSGrVJk66oUavaa1LvQVKjhmnSoEata6xJXQdJjZrXeZPUqActNanfIKlRJ7ptkhr1o5kmdRok\nNepKh01So9600aQeg6RGHeqqSWrUpwaa1F2Q1KhbnTRJjXpWe5P6CpIada75JqkRVTepoyCpEUPT\nTVIjinqb1EuQ1IgrTTZJjbiu0iZ1ESQ1YktjTVIj3lVjk9oPkhpxo2aapEZbVksPyFFdk44uLy+X\n3jAhNeJuF+dv//7JP558+8cvX75cess+1Oi64/O3z158vnpzPgzD5vHJZ598eHHyaOlRy7s4f3v2\n67/++Q+frVarpbfco+UTkhpxr6rPSWq05apGwzCs3pw/e/H5sntCVHROavaEtNls1uv16rEHlftd\n/Ps/x1/95unpaf5PyOuOjo5Wj09WH7y39JAUp797s/XN0189XmRJoIvztxd/+2/4Tf61pQdMZbVa\nnZ6eLr2CmiQ/qDd69erV0hOCfGsYhmE7SL/4xk//tcCWSO8Pw/vpN3mzJySgO0+fDmdnX14+eTL4\nVVqVlv9DAvry6afD8+dffFajCjkhAW3ZbIbsF1PcRpAAiOCVHQARBAmACIIEQARBAiCCIAEQQZAA\niCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIg\nARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAE\nQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmA\nCIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggS\nABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQ\nJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACI\nIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiAB\nEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARB\nAiCCIAEQQZAAiCBIAEQQJAAiCBIAEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAI\nggRABEECIIIgARBBkACIIEgARBAkACIIEgARBAmACIIEQARBAiCCIAEQQZAAiCBIAEQQJAAiCBIA\nEQQJgAiCBEAEQQIggiABEEGQAIggSABEECQAIggSABEECYAIggRABEECIIIgARDhfwQFVkMYm57C\nAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "elem = [2 3 1];\n",
    "subplot(1,2,1); showmesh(node,elem);\n",
    "[node,elem] = bisect(node,elem,1);\n",
    "subplot(1,2,2); showmesh(node,elem); \n",
    "findnode(node,5,'noindex');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Coarsening Procedure\n",
    "\n",
    "Our coarsening algorithm consists of finding good nodes and removing\n",
    "them. The coarsening procedure is illustrated by the following figures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBxE2xZjGwAAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNToxNzo1NPOqVNEAAB2z\nSURBVHic7d0xj+PImcbx8p2BnUhicoZxEAElBxirDRq4YE9KLEeesIGZ1aaLAzbvQX+A6fkAfdv5\nJh1eC2ugL+vMmkSCo1OiwwY2IIAKDEdsRXa0F9QMRyuJZBVZRb4k/7/AWC+4b1N6q+tpUWTVL376\n6ScFAEDd/qnuEwAAQCkCCQAgBIEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAg\nAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAA\nACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQC\nCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBE\nIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIA\nQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAg\nAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAI\nBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIA\niEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgk\nAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCB\nQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAA\nEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEE\nABCBQAIAiEAgAQBE+GXdJwDYeffuXd2n0Alv376t8sfR1mpU3FZbrQ2kxWLxu9/9ru6zaL8//vGP\n0+m0sh/37t27m5ubi9koCHtuK2+Xu+0qcl45jvbr+UYp5anycBwOJwOHZdXBOVc2edHWQ9vlTinV\ngrYW0NpAev/+fRD2Lmaj6fXEbeX7V/N49+y2chzt718/xNF+OA4v7146HOLbZfT45ikY9L/5w8xV\nTe3x6mm7it6/f19lICmlLmajONpPryduJ4KFWgZhz21l3dbhOPRXeTgZuB3hSWWHNU3Q1sRCLZVS\n7WirrTZ/h3QxGymlFrdL55Wn15M42ruqrMfK5XcvlVKXdy8fr57iaO+k8nYZLW5XurJbj1dP6uM7\nXLEg7Ll9lzxVTto6nAz8VS5frZrKuWirPzW21VabA0l9/CvDeSYFg56rTDoYK6Fy+suj08j5ByP1\nMY0u79znnCH5k9dRW71WdsVfZUNdbqs/tbfVSssDSfnLpNBBJp0dK06GeIvTSJM8eaVNAf4qlydk\n2upmW/0R0lZz7Q8kJTWTMsZKySHe+jTSZE5e2VOAv8pliJq2utZWf0S11VAnAknJy6TcsVJ4iHck\njTRpk5fJFOCvcjECp63utNUfgW010ZVAUpIyyXCsFBjinUojTc7kZT4F+KtsS+y01YW2+iO2rbk6\nFEhKRiZZjRWrId7BNNIkTF62U4C/yuaET1vtbqs/wtuarVuBpOrOpAJjxXCIdzaNtHonr2JTgL/K\nJhoxbbW1rf40oq0ZOhdIqr5MKjxWcod4x9NIq2vyKjMF+KucrUHTVvva6k+D2pqmi4Gk6sikkmMl\nY4iTRonqJ6/yU4C/ymkaN221qa3+NK6tZ3U0kFS1meRkrJwd4qTRkSonL1dTgL/Kpxo6bbWjrf40\ntK2nuhtIqqpMcjhWjoY4aXRWNZOX2ynAX+VDjZ62mt7W8tXSNLqtRzodSKqSTHI7VpIhvp5vSKM0\nFUxezqcAf5W1FkxbjW6rp0xqQVsPdT2QlM9MupiNFrer4Th0O1aCsHfx9ejx6slHZrQgjTSvE8F2\nFSmlgrDvvHK8ew4GfeeVWzNtNbGt/s65NW1NEEhK+cmkONo/vnnS/+C28nYZrR82SinnQ7w1aaT5\nmAjiZN+Hceij8uV3Lz2dc2umrWa1VfN3zq1pq0YgfeA2k5KxEoS96fXYYSZ92FHi7qXzId6yNNKa\nsvXA4eTir3L5akI0pa2HaKsJAukTV5l0NFaCsO9qr4qjuxgcDvFWppEmf+uBKiu3hvw3/xRtzUUg\n/Uz5TDo7VpzsVXH2njonQ7zFaaRJ3nqg+sqtIfnNT0NbsxFIx8pkUsZYKZlJGXd4lxzirU8jTebW\nA3VVbg2Zb3422pqBQDqjWCbljpXCmZT7vFHhId6RNNKkbT1Qb+XWkPbmm6CtaQik82wzyXCsFMgk\nw6dfCwzxTqWRJmfrAQmVW0POm2+Otp5FIKUyzySrsWKVSVZrMVgN8Q6mkSZh6wE5lVtDwptvi7ae\nIpCymGRSgbFimEkFVgYyHOKdTSOtiTtK1LtXRSM0cUcJ2nqEQMqRnUmFx0puJhVepy53iHc8jbQm\n7ihR114VDdLEHSVo6yECKV9aJpUcKxmZVHLV1IwhTholmrijRPV7VTROE3eUoK0JAsnIaSY5GStn\nM8nJGt5nhzhpdKSJO0pUuVdFQzVxRwnaqhFIpg4zKd653FHiMJMc7ijxaYjv9oo0StHEHSWq2aui\n0Zq4owRtVUr9su4T8CiO9ttl5LDgcBwublfbVbRdRdPrsVLKpP52FQVRL6fyJPywXuqbp8vvXuaW\n1Rlj8tP1uuBKqaEKp9djt29IHO3V5w7rGdkudwvlfmn29XyjlBqOQ93f3HNQSpmchm1lc/4qK/0C\nf++wntFPbFxbh+Ow9nO2Un1bbbU5kNbzjfO/YuLds/6H7XKnh2/O8dH+8eppOM7/M0dXjqP94nZl\ndCbGRyb1rY43sV1Flzdfua1p8kODsBeEORlfjHlZ29PwdMI+KsfRfruKqg6kBrZVyDm3SZsDKQh7\n8e756k/fuiqoP0cHYS8Y9IeTgb6Il+3uy++H41Bfl8s4TF+pi6P9xWyUe3ByJiZX9h6vnoYq1Lvs\n6DXCc/8TQ/ev5iZB69zFbBRH++n1xNVriT9uPRCEPcPKC7U0ObhAZX/nbFvZSTVztNVrZc355znn\n2vwd0sVsNL2e3H35vZNqyVXdYNDX1+sMV1u4mI2yDz783sjVuuCavlg3vR7rNHJ44fv+1Xx6PR5O\nBk6qWXF7Ef/gYr3d1gO5Bxeu7O+crSqXr2ZFyFvUyrY2SJsDSbnLpNPvGM3XcQjC/sXsi7SDT3eU\ncJVJR3cxOBziH9Ootq9bXb2W07a6mrxKVvZ3zlaVKybkLWpZW5ul5YGkXGRS2i+qTSb1zmZS2o4S\n5TPp7D11ToZ47WmklX8taW0tP3k5qezvnK0qV0zIW9SatjZO+wNJlcuk7F/UMpmUvaNEmUzKuMO7\n5BAXkkZamdeS3dYyk5fDyv7O2apyxYS8RS1oaxN1IpBU0Uwy+UUtlkkmO0oUy6Tc540KD3FRaaQV\ney0mbS02eW2XEXtVlCfkLWp0WxuqK4Gk7DPJ/BfVNpPiaP/45slkRwnbTDJ8+rXAEBeYRprtazFv\nq+3kNb0e37+eT68n7FVRnpC3qKFtba4OBZKyySTbX1SbvSqe9YNvhgFmnklWazFYDXGxaaSZvxbb\ntlpVfnzzpJRaPxg9/SbknGWmkSbkLXp883R597JBbW20bgWSMsukYr+oJpmkr9Rd3r3Uz3k4zKQC\nKwMZDnHhaaSZvJZibbWqrJRyeIdxNecsubNC3qKL2agpbW26zgWSysukMr+o2Zl0+L1R2n13Z+Vm\nUuF16nKHeCPSSMt+LWXaalXZ4VMvlZ2zZELeoka0tQW6GEgqPZPK/6KmZdLpXQyuMqnkqqkZQ7xB\naaSlvZbybbWq7GTyqvicJRPyFglvazt0NJDUuUxy9Yt6mklp99SVzyQna3ifHeKNSyPt9LW4aqtV\n5ZKTVy3nLJmQt0hsW1uju4Gkfp5Jbn9RDzNpu8q6w7tMJjlJo6Ty4RBvaBpph6/FbVutKheevGo8\nZ8mEvEUC21q+mhxtXlzVhF5o7u7L7+NofzEbxbu9vgUuQ7x73q4ivQFEhiDsL26XOjym15PDsnoP\niKMfpP9v7k9XSgVhTy9kvlXRUeWzzv64s4aTwf3rB71YXxPnrMSnJ0hWUcU7Shz9YRHY7FMQhL37\n1w9KKbeZ8Wny2j03NI20Gtt69mAhbb10umhyvboeSEqpw4WNDbcLWtyudJJlCwZ9/ZfRUdk42m9X\nURztT4fR49WTSeX449Pjhies18k3OTKO9sGgH4R9k4OFi3fPgc3uAP52lDCvrEeFj/dfr/jegs7S\n1kT7MqnrgZR8jo53+8Xt0mSvivtX82DQV3nXyvT3RkHYO91+Yj3fBGFvu4qSYaQ/8cTRXq8jnr39\nhA6t9XxjuFfFdhnF0d7kyt79q/k3P8yCsN/0IZ60Vb+WyrYeWNyuDtuRfXBGZbfv/9G70dzO1tVW\nVwf7ePNblkmd/g7p8Kqu1ToO0+uxHl5pBxztKKHOXcn55oev718/xNF+Pd+sHzb64Nzvkx4/7Cgx\nMXw+yVzyvVHTL08fttXh3bo/r5y/QYC/g60Ufjekoa1pGt3WI90NpNPvGC0zaZKWSaf31GVnUpJG\nKu8eh8eTHSVcZdLRXQzNHeKnbXU1eVlV9newFX+VK0ZbC7/AZuloIKXd8VI+k9Lu8D6bSdtVFAz6\n8e75cBilZdLZe+qcZNLZe+qaOMTT2lp+8rKq7O9gK/4qV4y2mmhcW8/qYiBl339ZJpOy1/A+yqTk\nSl1y7S458jSTMu7wLplJGXd4N2uIZ7e1zORlVdnfwVb8Va4YbTXXoLam6VwgZY8VrVgm5e4ooT5m\n0vphs13ukit1ybW7wyMPMyn3eaPCmZT7vFFThrhJW4tNXlZbDyilzA9mU4NcQtrq72Ay6Ui37rIz\nGd9a8nySyX130+vJ49XT/et5EPaOYkz/xNPTUEodHaz/5el/ro/MPQ2dSYvbpX7sKfeclfHTr/Jv\n4zFvq9VrCT5uPaCU0v97ZHG7SjkZo4PvX0eXd47PWXl7N6rnu62Gb76/gz29+cLbmq1DgWQ+vjXz\nTNL3VU+vx/HJ3dV3X35/9Jzper55vHoKwp7+HJacWxJF3/zwdTKM9DOAwaBvEjNWmWS1FoPkIW7b\nVvPXEh9sPXB0sP4kenqHdxzth+PQ/OD1w0bfEOzwnD29GxWrpq25b76/g63O2ZbYtubqyiU72/Gt\nmVy7S67Upd13N70eL25XyXIJenbTaaenKp1G3/wwU0odfp+kS13MRsPJQDndq6LAykAyLwUUa6vJ\na0kqm2w9kBysDLafsDq42Dk3vbPS2ur8YKtzLkZgW010IpCKjW8tO5OOvjc6m0lB2L+8e7m4XT1e\nPR3e4Z3c46DTKFnlXmfS4fdGQdh3uFdF4XXqpA3xMm3Nfi1Hles6uMw5WxHVWdrqiqi2Gmp/IJUZ\n31paJp29iyElk3rDyWA93+iFGA4qf6FXITpa5X44DtfzzeFFnrR7wc/KyKSSq6bKGeLl25r2Ws5W\nrv7g8udsRUhnaavt662rsictD6Ty41s7zaSMe+pOM2k932yXu6s/fZtcu1MHV+qCsHe6o8TVn77N\nvRc8w9lMcrKGt4Qh7qqtp68lo3KVB7s6Zyu1d5a2kkn/fHNzU/c5ePH+/fv1X/538V9LVyvs/nr0\nqxf9F//9n//zovdC/UKt5/+XcYf3cBL+dfO39XyjFwX/8ekv3/xh9qL/2XASPr1d/H3/j7/v//H4\n5klfqRtOPqwuvF3t9KC5vHv5ov/Zb17+2/3rhxe9F0HY//XoV0qpF/3PgrD/183ffnz6c7zb/8e3\n/55xwi/6n/169Ksfn/7849NflFL685mTt0K/EB2cv/mX0XQ6LV/TkPO2Hr6Wv+//kV059+DDtewK\nH6z/2dU5W6mrs8Lb6uRgq7bmHmylxl9YW22+y06vD2+yRr25YNDXBafX49xPKnpJ78Xt6vDgIGWV\n+2T7icOD9bW7ONrH0XNSNo72+r8yfOpIrwTh/q0Ie+v55vLmK1cFDXlpq/3WA2kHn91+wupg800N\nzM/ZSi2dFd7W8gdbtdXwYCt1/cJaaXMgKaXqveXxYjbSC34f/kv9FdHpfhD6s9HxB9bP1ZkB9LlS\nnxufhNXBVj5XlzdfvX371k/1LP7aalI5rYMZptPp6Z+l05vfHx93tt1prA62UlNnnbfVqlOu2npG\nt9tqpc2BdDEbxXmLw1vZLqPtcqefMFjPNybPzOo/do6eT1qHmyDsxdE+uZ3h8eOOEsKHiwTO2xrb\n7CagpW1qkLb9xHQ6pbPZnLdVM9x+wupg2upPm29qcPtt3uFdDOXXYNXnpu9xyF0ZCIfctvXgS2m7\nDQJyDz6sXP48W0/It/q0tV5tDiTlbpSf3lPnKpP0MjOkkRVXbT29RcrV5OXqhrFOkZ9JtNW3lgeS\ncjHK0+7wLp9Ji9ulvlJhuBM5EuXbmja5lJ+8mLYKk5xJtLUC7Q8kVW6UZ6/hXSaTkit1ybU723Pr\nuDJtzZ5cykxeTFslycwk2lqNTgSSKjrKTXaUKJBJi9ulvmf048pAH75PMj8xaMXaajK5FJu8lMH2\nE8glLZNMdpSAE22+y+5IMrwMV8A1SSPNdq8K/WjRVkWGMYYMtm01/1PXqnLwYeuBKD63/cSZO7yR\nybatnioHNjtKoLwOBZKyGYvmaaSZZ9Lj1dNwHAazXnyyV8XNv94a/jgcMm+r7YUXq8rZe1XAloRM\nii13lEBJXblklzD5zG6bRprJtbvke6O0vSpQjElbi30NYFXZausB5Kr32h1trV7nAknljcViaaRl\nZ9LR80ZkklvZbS3zpbRV5WatZSlfXZlEW2vRxUBS6cOrTBppaZl09ulXMsmttLaWv0XKqjKTl1vV\nZxJtrUtHA0mdG17l00g7zaSMtRjIJLdO2ypkUwOUUWUm0dYadTeQ1M+Hl6s00g4zKXdlIDLJrcO2\nun18xKoyk5db1WQSba1Xt+6yO/XpUYNVdDEbmaTCdhWt5x8Wn8+lh/hFmF/ZsCBMfJo1ds9uHx/5\nuODTg+GOHkET1vxvCt/33dHW2nU9kJRSequhIOwFYS8I+wbH75VSJnNcMqZNDk520IIr8e45GPRN\nelqAHjAmhw3HPL/ijL9MSurT1hp1PZCSK3VxtH+8erq8+yJ3OK4fNsm2ERn080ZbFV3MRovbpckz\nswSSK8mFlyDsu528jiobbmrg5EdD85FJtFWITn+HdPi9kdtLw4ffG13MvjBfWwjlHX4N4G2vCseV\nYYW2tlV3A+n0LgZXY/H0Lgar9e5QRskdJWqpjAJoayt1NJDS7qkrPxbT7qkjkyrgZEeJiiujMNra\nPl0MpOw7vMuMxew7vMkkrxzuKFFZZZREW1umc4Fk8rxRsbGY+7yRIpO8cb6jRAWV4QRtbZNuBZL5\n06+2Y9EkjTQyybkCO0oYttVfZThEW1ujQ4FkuxaD+Vg0TyONTHKo8I4SuW31VxnO0dZ26EogFVsZ\nyGQs2qaRRiY5UcGOEs4rwxPa2gKdCKQy69Rlj8ViaaSRSSVVtqOEw8rwirY2XfsDqfyqqWljsUwa\naWRSYRXvKOGkMipAWxut5YHkag3v07FYPo00MqmAWnaUKFkZlaGtzdXmtey2y912uXO1o8ThAtKL\n25X+v04qX8xGSikyyZCnHSX02mW+Kn9evh4s0NaGanUgraLhOLx/NXdYM9496x0lhsqosl6z1eEJ\nYD3fmGwQYEXvJqCU8lSZfQpybZc758uV0tbGaXMgXcxG+sOHQ4vblf60fvH1KBjkrzT8+ObZ+Ehy\ny8hwHA4nA0/FfexogNo5byvbT3jS5kAKwp7DndmUUvpz+lCFF1+P1g8bw9Xvh+Mw97D7V/PL717e\nv3b5Ya6thpPB9HrisKC+BqjbFEd7k60HbCs7qdZuDWqrxvYTPrT8pgaHDu9iCAYuv9u8fzWfXo/d\nZicMHXwjNfC2qYGvz3NI46+t8IpAMnJ6T52rUU4a1ajKvSpQGXaUaC4CKV/aHd7lRzlpVKPq96pA\nBdhRotEIpBzZzxuVGeWkUY3q2qsCXrGjRNMRSFlMnn4tNspJoxrVu1cFPGFHiRYgkFKZr8VgO8pJ\noxpJ2KsCzrGjRDsQSOfZrgxkPspJoxrJ2asCDrGjRGsQSGcUW6fOZJSTRjWStlcFnGBHiTYhkI6V\nWTU1e5STRjWSuVcFSmJHiZYhkH6m/BreaaOcNKqR5L0qUBg7SrQPgfSJqx0lTkc5aVQj+XtVoAB2\nlGglAukDV2mkJaNcKXX/+oE0qounvSr0iu+kUV38tbV8NZRBICnlOo00Pcr1qo7MWbXwkRnJ5EUa\n1cVrW8mkehFIXtIoqXx593L9sGGUV8/rJ5h49xwM+kHYd14Z2fy1lUySoOuB5C+N9PdGF7MRo7x6\n/qatpDJtrZ7vy6RkUu06HUi+00j/2jDKK1ZBGg0nIW2tWDVf2tHWenU3kKpJI41RXplq0kj/G9pa\nmSpvIaGtNepoIFWZRhqjvAJVppFGWytQ/Q2NtLUuXQyk6tNIY5R7VX0aabTVq7pur6ettehcINWV\nRhqj3JO60kijrZ7U+7AXba1etwKp3jTSGOXO1ZtGGm11TsKjx7S1Yh0KJAlppDHKHZKQRhptdUhC\nGmm0tUpdCSQ5aaQxyp2Qk0YabXVCThpptLUynQgkaWmkMcpLkpZGGm0tSVoaabS1Gu0PJJlppDHK\nC5OZRhptLUxmGmm0tQItDyTJaaQxyguQnEYabS1AchpptNW3NgfSer5RnlZNffPkcEeJw70qkEt+\nGmlMXlbkp5FGW736Zd0n4Mtvf/vbm5ubdbTRseTc/eu585rT6dR5zfZZzzfDcbhdRdtV5LDsdrnb\nriLnlYOwt55vLm++clWwrfy1VSm1UEuHNRVt9eYXP/30U93nAFh49+5d3adg7e3bt3WfgnS0FYpA\nAgAI0ebvkAAADUIgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAA\nACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQC\nCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBE\nIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIA\nQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAg\nAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAI\nBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIA\niEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgk\nAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCB\nQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAA\nEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEE\nABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQ\nSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAg\nAoEEABCBQAIAiEAgAQBEIJAAACL8P8HXhoN3STalAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "load Lshapemesh\n",
    "figure\n",
    "for k = 1:3\n",
    "    subplot(1,3,k); showmesh(node,elem);\n",
    "    [node,elem] = coarsen(node,elem,'all');\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 1 -- Find Good Nodes\n",
    "We use the following characterization to find a good node p\n",
    "\n",
    "- The valence of p is 4 (interior node) or 2 (boundary node).\n",
    "- It is the newest vertex of all triangles in its nodal star."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBxsJiREDdwAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNToyNzowObfMpvkAACAA\nSURBVHic7d09bxtZloDh0qJjquKFCuDGWmC5kUElZjbKWgN76bDdWDg2e/kD3P0DCMuxg7FDE26M\nnCkbOiHhjEEz56CEiSn+AW1wbDatD7I+7se5t94nc2MkXssevi7WrXMPbm5uEgAAfPs33wsAACBJ\nCBIAQAmCBABQgSABAFQgSAAAFQgSAEAFggQAUIEgAQBUIEgAABUIEgBABYIEAFCBIAEAVCBIAAAV\nCBIAQAWCBABQgSABAFQgSAAAFQgSAEAFggQAUIEgAQBUIEgAABUIEgBABYIEAFCBIAEAVCBIAAAV\nCBIAQAWCBABQgSABAFQgSAAAFX7wvYD4TSaTz58/+14FAG8eP37c6/V8ryIABzc3N77XELmDg4M0\na3X6x74X4scqX8/Hi8b+BOS3nyRJp3+cZi3fy/FgOb1azvLG/vaTJJmMZkmS8E5bBFdILsh7cW94\n4nshHiyn+XKWt7tZ0ryfwCpfv3v6od3N0qy1yte94UnT3pSX03w5vWp3s07/uH2S+V6OBxeDy07/\nWP5Rgr24h+SCvBFPRlPfC/EjPTrsDU9W+bpRPwGp0dnr0/bJUZq1zs5PLwaXq3zte13uLKf5ZDR7\n/nvf90K8uRhcJklydn7qeyHBIEiONL1JWatRTdqq0dfLgqY1iRpRowoIkjs0qSFNulsj0ZwmUSNq\nVA1BcoomRd+kh2okmtAkakSNKiNIrtGkiJu0u0Yi7iZRI2pUB0HygCZF2aQiNRKxNokaUaOaCJIf\nNCmyJhWvkYivSdSIGtVHkLyhSdE0qWyNRExNokbUyAiC5BNNiqBJ1Wok4mgSNaJGphAkz2hS0E2q\nUyMRepOoETUyiCD5R5MCbVL9Golwm0SNqJFZBEkFmhRck0zVSITYJGpEjYwjSFrQpICaZLZGIqwm\nUSNqZANBUoQmBdEkGzUSoTSJGlEjSwiSLjRJeZPs1UjobxI1okb2ECR1aJLaJtmukdDcJGpEjawi\nSBrRJIVNclMjobNJ1Iga2UaQlKJJqprkskZCW5OoETVygCDpRZOUNMl9jYSeJlEjauQGQVKNJnlv\nkq8aCQ1NokbUyBmCpB1N8tgkvzUSfptEjaiRSwQpADTJS5M01Ej4ahI1okaOEaQw0CTHTdJTI+G+\nSdSIGrlHkIJBk5w1SVuNhMsmUSNq5AVBCglNctAknTUSbppEjaiRLwQpMDTJapM010jYbhI1okYe\nEaTw0CRLTdJfI2GvSdSIGvlFkIJEk4w3KZQaCRtNokbUyDuCFCqaZLBJYdVImG0SNaJGGhCkgNEk\nI00KsUbCVJOoETVSgiCFjSbVbFK4NRL1m0SNqJEeBCl4NKlyk0KvkajTJGpEjVQhSDGgSRWaFEeN\nRLUmUSNqpA1BigRNKtWkmGokyjaJGlEjhQhSPGhSwSbFVyNRvEnUiBrpRJCiQpP2NinWGokiTaJG\n1EgtghQbmrSjSXHXSOxuEjWiRpoRpAjRpHub1IQaiYeaRI2okXIEKU406VaTmlMjcbdJ1Iga6UeQ\nokWTNk1qWo3EdpOoETUKAkGKGU2SJp0/etu0GolNk6hRQo1CQJAi1/AmJUkyHy/SrLWc5b4X4scq\nv17O8tXVtbOzz7WhRgEhSPFrbJPkk7pO/zg9OnRz9rk28kldp3/c6R+7OftcG2oUFoLUCA1s0vZ9\nIzdnn2uzfd8ozQ4dnH2uDTUKDkFqikY16e4uhqY16e4uBttnn2tDjUJEkBqkIU16aE9dc5r00J66\n5jSJGgWKIDVL9E3avcO7CU3avcO7CU2iRuEiSI0TcZOKPG8Ud5OKPG8Ud5OoUdAIUhNF2aTiT7/G\n2qTiT7/G2iRqFDqC1FCRNansLIb4mlR2FkN8TaJGESBIzRVNk6pNBoqpSdUmA8XUJGoUB4LUaBE0\nqc6cujiaVGdOXRxNokbRIEhNF3ST6k9NDb1J9aemht4kahQTgoRQm2Rqhne4TTI1wzvcJlGjyBAk\nJEmATTJ7okSITTJ7okSITaJG8SFI+CqgJtk43yisJtk43yisJlGjKBEk/CmIJtk7bS+UJtk7bS+U\nJlGjWBEkfEd5k2yf/aq/SbbPftXfJGoUMYKE29Q2yc1J5Jqb5OYkcs1NokZxI0i4h8ImuamR0Nkk\nNzUSOptEjaJHkHA/VU1yWSOhrUkuayS0NYkaNQFBwoOUNMl9jYSeJrmvkdDTJGrUEAQJu3hvkq8a\nCQ1N8lUjoaFJ1Kg5CBL28NgkvzUSfpvkt0bCb5OoUaMQJOznpUkaaiR8NUlDjYSvJlGjpiFIKMRx\nk/TUSLhvkp4aCfdNokYNRJBQlLMmaauRcNkkbTUSLptEjZqJIKEEB03SWSPhpkk6ayTcNIkaNRZB\nQjlWm6S5RsJ2kzTXSNhuEjVqMoKE0iw1SX+NhL0m6a+RsNckatRwBAlVGG9SKDUSNpoUSo2EjSZR\nIxAkVGSwSWHVSJhtUlg1EmabRI2QECTUYaRJIdZImGpSiDUSpppEjSAIEmqp2aRwayTqNyncGon6\nTaJG2CBIqKtyk0KvkajTpNBrJOo0iRphG0GCARWaFEeNRLUmxVEjUa1J1Ai3ECSYUapJMdVIlG1S\nTDUSZZtEjXAXQYIxBZsUX41E8SbFVyNRvEnUCPciSDBpb5NirZEo0qRYaySKNIka4SEECYbtaFLc\nNRK7mxR3jcTuJlEj7ECQYN69TWpCjcRDTfJSozRfn+XrTn7t9EUfaBI1wm4ECVbcalJzaiTuNslL\njdrT/OWjt69m+avR7OWjt6nD4/XuNokaYS+CBFu+NWmWJEmjaiQ2TZqPF8vplZdP6n56Ov5zPfn6\nx8Gly1f/s0lX1/LXgBpht4Obmxvfa4jcwcFBu9ugN+JblrM8SZI0a6VHh14WsLq6TpLE46vLJYL7\nvwOdq+u/f39JtMpa/+385+DxJ6DHcpbzTlvED74X0AjLWf78Y8z3sR+yulrL+1F6dNgbdr2sYTnL\n5+NFkiReFrCc5XJx4H4B7SRJtq6QkiRJ83XvtetrlMloJkHqPDtOj1qOX12Di1+cXpgGjSC5cHZ+\nevHL5eDLC98LcWqVry9+uez0j5fTq/bJ0XKWy4d4rpdxtU6PDr0sYDnNl9Or3rAr78geFtDN2rN8\n88tP564/Nb0YXKZZq51knWfH8w+Ls/PTNGtWk949GZ+9Pn33/b8M8BDuIbnQ6R/3hifnj976Xog7\nf+5i6GaJk7PPd0izVqf/n44XsL2Lwc3Z53d9Oj+dDLvzrLXKWpNhd94/dvnq27sY0iMXZ59r8+7J\nuDfsNurWaU0EyZFGNenePXWNatLdPXVemrTKWp+HJz93s5+HJ5/dXpzd3VNn++xzbahRBQTJnYY0\naccO74Y06aEd3r6uk9x7aId3c5pEjaohSE5F36S9zxtF36Tdzxs1oUm7nzdqQpOoUWUEybWIm1Tw\n6deIm1Tk6de4m1Tk6de4m0SN6iBIHkTZpFKzGKJsUvFZDLE2qfgshlibRI1qIkh+RNakCpOBImtS\n2clA8TWp7GSg+JpEjeojSN5E06TKc+qiaVK1OXUxNananLqYmkSNjCBIPkXQpJpTUyNoUp2pqXE0\nqc7U1DiaRI1MIUieBd0kIzO8g25S/RneoTep/gzv0JtEjQwiSP4F2iSDJ0oE2iRTJ0qE2yRTJ0qE\n2yRqZBZBUiG4Jhk/3yi4Jpk93yjEJpk93yjEJlEj4wiSFgE1ydJpewE1ycZpe2E1ycZpe2E1iRrZ\nQJAUCaJJVs9+DaJJ9s5+DaVJ9s5+DaVJ1MgSgqSL8iY5OIlceZNsn0Suv0m2TyLX3yRqZA9BUkdt\nkxzUSKhtku0abRagtkm2ayQ0N4kaWUWQNFLYJGc1Egqb5KZGmwUobJKbGgmdTaJGthEkpVQ1yXGN\nhKomuazRZgGqmuSyRkJbk6iRAwRJLyVN8lIjoaRJF4NLxzXaLEBJk9zXSOhpEjVygyCp5r1JHmsk\nvDcpzQ7n40X75MjXArw3yVeNhIYmUSNnCJJ2HpvkvUbCY5OW03z+YZFmLV8LSHw3yW+NhN8mUSOX\nCFIAvDRJSY2ElyZt3zdyc/b5Q3w1SUONhK8mUSPHCFIYHDdJVY2E4ybd2sXg4Ozz3dw3SU+NhPsm\nUSP3CFIwnDVJYY2Esybdu6euUU3SViPhsknUyAuCFBIHTVJbI+GgSTt2eDekSTprJNw0iRr5QpAC\nY7VJymskrDZp7/NG0TdJc42E7SZRI48IUngsNSmIGglLTSr49GvETdJfI2GvSdTIL4IUJONNCqhG\nwniTSs1iiLJJodRI2GgSNfKOIIXKYJOCq5Ew2KQKk4Eia1JYNRJmm0SNNCBIATPSpEBrJIw0qfKc\numiaFGKNhKkmUSMlCFLYajYp6BqJmk2qOTU1giaFWyNRv0nUSA+CFLzKTYqgRqJyk4zM8A66SaHX\nSNRpEjVShSDFoEKToqmRqNAkgydKBNqkOGokqjWJGmlDkCJRqkmR1UiUapLx842Ca1JMNRJlm0SN\nFCJI8SjYpChrJAo2ydJpewE1Kb4aieJNokY6EaSo7G1SxDUSe5tk9ezXIJoUa41EkSZRI7UIUmx2\nNCn6GokdTXJwErnyJsVdI7G7SdRIM4IUoXub1JAaiXub5KBGQm2TvNQozdcvr657s9zpiz7QJGqk\nHEGK060mNapG4laTnNVIKGySlxp1xouXj94O8vVPo9nLR29Th8fr3W0SNdKPIEVr06QG1khIk+Yf\nFsuZ0xqJTZPm44XL191egDRpOcslS+6vjX4cXD70SxcL2GoSNQrCwc3Nje81RO7g4ODXfw19vfpk\nNJ2MZmnW8rWAVb72+OqygCRJ6qyhzm+h/qvXtLk+cL+Gs3z9t1uLyVr/4XgR334CZ+ennf6x8xf/\n6td/H/FOW8QPvhcAi1b5ej5edPrH8/Fi8OWFpzVcv3s6bnczLzfSl7P8YnCZZq12N5MLpgrePf2w\nytfPP/bT7LDs187Hf0xGsyRJ6iygjsloKpdo7hfQya+Tp+Pt/zI/Onzu/K/BxeCy8+x4/mHR7mZ+\n/22EvQhStLY/qWufZO+efvDSpFV+3e5m7ZOjyWjquEnLaT7/sDg7P11O8zRrzcd/VH5HHnx58e7p\nh+cfn5V6R5uPF8vpVW/YlWusOguoRu4bdfrHadZa5WvHC1hmrWU3a29tZ/jnsOs4CZtP6trd7GJw\neXZ+SpM04x5SnG7dN3Jw9vluveFJmrUuHN5CuLWLoeYM1jRrPf/4TC6VCn7JfLyYf1jIArzscdje\nxZBmhw7OPr/r0/npZNidZ61JN3v/sb90ewtn+76Rm7PPURNBitC9uxga1aR799S5bNJ2jTZf7rJJ\nd/fU2T77/F6rrPV5ePLXo8P3w67HGgmapB9Bis2OPXUNadKOHd5umnS3Rpsvd9Okh3Z4e2mSFw/t\nqaNJyhGkqOzd4R19k/Y+b2S7SQ/VaPPltpu0+3mjJjRp9w5vmqQZQYpHweeNIm5Swadf7TVpd402\nX26vSUWefo27SUWeN6JJahGkSJR6+jXKJpWaxWCqSck0/3Fw+dOT8ePRtEiNNl9uo0nFZzHE2qTi\nT7/SJJ0IUgwqzGKIrEkVJgPVb9Kvw5NXT8ed8aI9y3uj2evBZfEFGG9S2clA8TWp7CwGmqQQQQpe\n5clA0TSp8py6mk16PP1uYGg7SdrTEiNEDTap2py6mJpUbTIQTdKGIIWt5py6CJpUc2pqnSbdHRWa\nXpV7XzPSpDpTU+NoUp05dTRJFYIUMCNTU4NukpEZ3pWbtLrzzP+yW/oPomaT6s/wDr1J9aem0iQ9\nCFKoDM7wDrRJBk+UqNakN1lrufXL8+9/WVzlJpk6USLcJpma4U2TlCBIQTJ+okRwTTJ+vlHZJs3H\ni8n06tOXFzId582XF/8sOVtoW4UmmT3fKMQmmT1RgiZpQJDCY+l8o4CaZOm0veJN2uzwXmWtef/4\n8/BklbUqzLvbVqpJNk7bC6tJNs43okneEaTAWD1tL4gmWT37tUiTdjxv5KZJ9s5+DaVJ9k7bo0l+\nEaSQODj7VXmTHJxEvrtJe59+td0k2yeR62+S7bNfaZJHBCkYzk4iV9skBzXaLCC5LwkFZzHYa5Lt\nGm0WoLZJbk4ip0m+EKQwOKuRUNgkZzXaLCD5PgnFJwMldprkpkabBShskpsaCZrkBUEKgOMaCVVN\nclyjzQKSb0koVSNhtkkua7RZgKomuayRoEnucYS5dl5qJDr94yRJzh+99XL2eZIkveHJZDR992S8\nurruDU/m40WFb7Kc5qt8Xe1r0+xwMprKm3K1BXT6x5PRrNM/rraA5Ns56PLmWO07VCZNmoym8tt3\n/Orb3NdIbJrE2eduECTVPNZIeG9Su5vJW/kqv672HVb5ejnL06xV8w2l2gLkH9fyW6i2gFW+XuXr\n3km3wtfWp6FJvmokaJJLBEkv7zUSHpskn9QNvryQy5RqlwjzTHq27vSPy/4k5+NFenSYfPshlH1H\nXuXr80dv290szVrLWf58+KzsO9rF4HJzdeUrCX6b5LdGgiY5wz0kpZTUSHi5n7R936jmDFZ5Q5mM\nZssy07i37xtVmC0kNXr+sd8+Oap2P2lz38jZ2ecP8XU/SUONBPeT3CBIGqmqkXDcpLu7GBw36e4u\nhlJN2qpRtllAqSbd2sXQwCbpqZGgSQ4QJHUU1kg4a9JDe+qcNemhPXUFm3S3RpsFFGzSvXvqGtUk\nbTUSNMk2gqSL2hoJB03avcPbQZN27/De26SHarRZwN4m7djh3ZAm6ayRoElWESRFlNdIWG1SkeeN\nrDapyPNGO5q0u0abBexo0t7njaJvkuYaCZpkD0HSIogaCUtNKv70q8Empfm6Pc3l6PHiT79uNynN\n153xIs3XRWq0WcCmSdsLKPj0a8RN0l8jQZMsYdu3CgHVSBjfC152FoNsRJaduBVeTt5QVoPLl7Ov\n10mrrPXXo8OyC/jx0dtOvpYvf5Ovi9Ros4DnH58lTz+8/PaOtkySi/5xwd+ONGk+/iOmveCh1Eiw\nF9wGrpD8C65GwuB1UrXJQPWvk/42+/NTuzRf/63kd3idrzvfcpLm65dZq1dyAa+3/n3dTpJSC4js\nOimsGgmuk4wjSJ4FWiNhpEl15tTVaVL7zj2k9KrcLIb27LvvkObr9KrEG1PnzjChW99wr2iaFGKN\nBE0yiyD5FHSNRM0m1Z+aWrlJq+yw8ot+/Q5Htb7D6uj25zzLbum/BhE0KdwaCZpkEEHyJoIaicpN\nMjXDu1qTVlnrVgA+l7wX8unZ8fYvl0ky7x8/8L+9dwGH8+/vPXyu9Dch6CaFXiNBk0whSH5EUyNR\noUlmT5So0KT5ePFzkkyG3VXWWmWtybA7+LAoPltola8Hg8v3w+6ym62y1rx//POwW+od+d1o+ls3\n2yzgzfnpr6Op1bPP7anWpDhqJGiSEQTJg8hqJEo1ycb5RqWaJDu8e7/3Pw9P3nx58ebLi8/Dk+Kz\nhTY7vJfDk/e/9998efHp/LTUbCFZZ+f8dLOAVf/Y6tnntpVtUkw1EjSpPoLkWpQ1EgWbZO+0vYJN\neuh5o4KzhXY8b1SwSQ89b2Tv7HM3ijcpvhoJmlQTQXIq4hqJvU2yffbr3ibtfvp1b5P2Pv26t0m7\nn35tQpNirZGgSXUQJHeir5HY0SQ3J5HvaFKRWQw7mlRwFsOOJhWZxRB3k+KukaBJlREkRxpSI3Fv\nk9zUSNzbpOKTge5tUvHJQMkDTSo4GSiJt0lNqJGgSdUQJBcaVSNxq0kuayRuNal4jcStJpWq0WYB\nyVYSitdos4DImtScGgmaVMHBzc2N7zVE7uDgIM1ajarRxny8uBhcplkrLTMmzqDJaCrnf1dbwCpf\nXwwul7M8SZJSNbq7gHY3qzB2T/4pkyRJtZmBq3w9H/8xHy96w5NOmWekTJEgzceLxv5fYPNXiHfa\nIgiSdQcHB76XAMCnXq/3j3/8w/cqAsC0bxfkX8cNHAm8nOYXv1yu8nWatUzNBS9rPl5MRtNqFyji\n/NHb9Oiw2sdN8q/jwf/+36tXr6q9OtAc3ENyoZkfJct9o7PXp+1u5ubs84e0u1mdueBJkhR/Znbb\n16dffXxWBoSIILnQwNubt3YxODj7fDcHZ5/fUnYXAwCC5EijmnTvnrpGNYkaARUQJHca0qQdO7wb\n0iRqBFRDkJyKvkl7nzeKvknUCKiMILkWcZMKPv0acZOoEVAHQfIgyiaVmsUQTZOS0fSnJ+OXj97+\nOLicUyOgHoLkR2RNqjAZKIImDZ4dvxrN2rM8zded8eL1eEGNgDoIkjfRNKnynLrQm/T4+4/s2knS\nLvmgEoBtBMmnCJpUc2pq0E1K7/yppVeh/jkCGhAkz4JukpEZ3uE2aXlydPu/dBs3PBQwiCD5F2iT\nDJ4oEWiT3uTr+dZ8wt+62Ty/Nr00oEEIkgrBNcn4+UbBNelicLlMkk9fXrz/2H//sf/bv4arSvPu\nAGwQJC0CapKl0/YCatL280bLk2x5kiWV5t0B2EaQFAmiSVbPfg2iSTuefqVJQB0ESRflTXJwErny\nJu2dxUCTgMoIkjpqm+SgRkJtkwpOBqJJQDUESSOFTXJWI6GwSaXm1NEkoAKCpJSqJjmukVDVpApT\nU7eadGVriUBcfvC9ADxo06Sz89N063kXx7zUSMjh3+eP3g6+vHD/6kmS9IYnk9F0NV7P80Vv2J2M\npmW/Q5q15uNF8hcbqwNiQ5BU894kjzUS3pu0ytftbra6uk6zwwpfXu2rgGYiSNp5bJL3GgmPTZJP\n6p7/3p+MpstpXm2Y94rxDUAx3EMKgJf7SUpqJLzcT9q+b1RzLjiAIghSGBw3SVWNhOMm3d3FQJMA\n2whSMJw1SWGNhLMmPbSnjiYBVhGkkDhoktoaCQdN2r3DmyYB9hCkwFhtkvIaCatNKvK8EU0CLCFI\n4bHUpCBqJCw1qfjTrzQJsIEgBcl4kwKqkTDepLKzGGgSYBxBCpXBJgVXI2GwSRUmAyU0CTCNIAXM\nSJMCrZEw0qRqNRI0CTCIIIWtZpOCrpGo2aQ6NRI0CTCFIAWvcpMiqJGo3KT6NRI0CTCCIMWgQpOi\nqZGo0CRTNRI0CaiPIEWiVJMiq5Eo1SSzNRI0CaiJIMWjYJOirJEo2CQbNRI0CaiDIEVlb5MirpHY\n2yR7NRI0CaiMIMVmR5Oir5HY0STbNRI0CaiGIEXo3iZ5qVGar38cL/5+df24/OHfddzbJDc1EjQJ\nqIAgxelWk7zUqD3NXz56ezZedPJ1bzT76cnY5avfapLLGgmaBJTFEebR2jSp8+x4Mpo+//jM5YGz\nSZL813ix/cv2LE/Hi2U3c7aAdjfr9I8lCe1udnZ+6vgn0On/52Q0nYxmvV//4vJ1gUARpJilWavz\n7PhicJlmrXdPPzh+9Zd3b2KNpu8cL+Kb1dW1+58AgFIIUsyW03z+YTH48uJicHl2fppmLacv/2Sc\nzPLt/9D++Gzgdg0Xg8tl8nUNgy8vXL70xsTt/TMgXNxDitbmvpGzs89v+XR+utrKz2TYXTmvUZIk\nveFJu5u5OfscQB0EKU63djF4adIqa73/+Oz9sHuetd58efF5eOLspZM7uxgcnH0OoCaCFKF799T5\natKkm705OvRybXRrTx1NApQjSLHZscPb12d3ju3Y4U2TAM0IUlT2Pm8UfZP2Pm9EkwC1CFI8Cj79\nGnGTCj79SpMAnQhSJErNYoiySaVmMdAkQCGCFIMKk4Eia1KFyUA0CdCGIAWv8py6aJpUeU4dTQJU\nIUhhqzk1NYIm1ZyaSpMAPQhSwIzM8A66SUZmeNMkQAmCFCqDJ0oE2iSDJ0rQJEADghQk4+cbBdck\n4+cb0STAO4IUHkun7QXUJEun7dEkwC+CFBirZ78G0SSrZ7/SJMAjghQSByeRK2+Sg5PIaRLgC0EK\nhoMaCbVNclAjQZMALwhSGJzVSChskrMaCZoEuEeQAuC4RkJVkxzXSNAkwDGCpJ2XGgklTfJSI0GT\nAJcIkmoeayS8N8ljjQRNApwhSHp5r5Hw2CTvNRI0CXCDICmlpEbCS5OU1EjQJMABgqSRqhoJx01S\nVSNBkwDbCJI6CmsknDVJYY0ETQKsIki6qK2RcNAktTUSNAmwhyAporxGwmqTlNdI0CTAEoKkRRA1\nEpaaFESNBE0CbCBIKgRUI2G8SQHVSNAkwDiC5F9wNRIGmxRcjQRNAswiSJ4FWiNhpEmB1kjQJMAg\nguRT0DUSNZsUdI0ETQJMIUjeRFAjUblJEdRI0CTACILkRzQ1EhWaFE2NBE0C6iNIHkRWI1GqSZHV\nSNAkoCaC5FqUNRIFmxRljQRNAuogSE5FXCOxt0kR10jQJKAyguRO9DUSO5oUfY0ETQKqIUiONKRG\n4t4meanR49H01Wj6j1n+eDR1+bo0CaiAILnQqBqJP5t0tU481agzXvRGs06+TvN1bzT76cnY6avT\nJKCkH3wvoBEufrk8e326nOa+F+Ja59mxpKidZL1h1/FP4OX3V0XtWZ6OF/OjlrMFpEetdjebjGa9\nX//i7EWBcBEkF9Kjw8lo5nsVPq2urt3/BNK7N7E+LJaOFwGgMILkQqM+rNt2MbhsJ9nq6jo9Ojw7\nP00zd1cnSZIsn4zbs++uyXrO/yDePRm3u5njFwUCxT0k2CIf1vWGXamRg7PPb/l0frr6lsBV1nr/\n0UONesNu++TI8esCgSJIsOLWLgYHZ5/ftcpa7z8++23Y/a2bvfnyYnni9ErlW424PAKKIkgw7949\ndb6aNM8OL9x+VJhQI6ASggTDduzw9tIk96gRUA1Bgkl7nzeKvknUCKiMIMGYgk+/RtwkagTUQZBg\nRqlZDFE2iRoBNREkGFBhMlBkTaJGQH0ECXVVnlMXTZOoEWAEQUItNaemRtAkagSYQpBQnZEZ3kE3\niRoBBhEkVGTwRIlAm0SNALMIEqowfr5RcE2iRoBxBAmlWTptL6AmUSPABoKESNFC1wAABG9JREFU\ncqye/RpEk6gRYAlBQgkOTiJX3iRqBNhDkFCUgxoJtU2iRoBVBAmFOKuRUNgkagTYRpCwn+MaCVVN\nokaAAwQJe3ipkVDSJGoEuEGQsIvHGgnvTaJGgDMECQ/yXiPhsUnUCHCJIOF+SmokvDSJGgGOESTc\nQ1WNhOMmUSPAPYKE2xTWSDhrEjUCvCBI+I7aGgkHTaJGgC8ECX9SXiNhtUnUCPCIIOGrIGokLDWJ\nGgF+ESQkSVA1EsabRI0A7wgSwquRMNgkagRoQJCaLtAaCSNNokaAEgSp0YKukajZJGoE6EGQmiuC\nGonKTaJGgCoEqaGiqZGo0CRqBGhDkJooshqJUk2iRoBCBKlxoqyRKNgkagToRJCaJeIaib1NokaA\nWgSpQaKvkdjRJGoEaEaQmqIhNRL3NokaAcoRpEZoVI3EVpOuE2oEhIAgxa+BNRLSpPl4MR8vqBGg\n3w++FwC7GlsjkWat9OgwPTpMs0PfawGwB1dIMWt4jZJvn9S5PPscQGUEKVrUaHPfyNnZ5wDqIEhx\noka3djHQJEA/ghQhanTvnjqaBChHkGJDjXbs8KZJgGYEKSrUaO/zRjQJUIsgxYMaFXz6lSYBOhGk\nSFCjUrMYaBKgEEGKATWqMBmIJgHaEKTgUaPKc+poEqAKQQobNao5NZUmAXoQpIBRIyMzvGkSoARB\nChU1MniiBE0CNCBIQaJGxs83okmAdwQpPNTI0ml7NAnwiyAFhhpZPfuVJgEeEaSQUCMHJ5HTJMAX\nghQMauSgRoImAV4QpDBQI2c1EjQJcI8gBYAaOa6RoEmAYwRJO2rkpUaCJgEuESTVqJHHGgmaBDhD\nkPSiRt5rJGgS4AZBUooaKamRoEmAAwRJI2qkqkaCJgG2ESR1qJHCGgmaBFhFkHShRmprJGgSYA9B\nUoQaKa+RoEmAJQRJC2oURI0ETQJsIEgqUKOAaiRoEmAcQfKPGgVXI0GTALMIkmfUKNAaCZoEGESQ\nfKJGQddI0CTAFILkDTWKoEaCJgFGECQ/qFE0NRI0CaiPIHlAjSKrkaBJQE0EyTVqFGWNBE0C6iBI\nTlGjiGskaBJQGUFyhxpFXyNBk4BqCJIj1KghNRI0CaiAILlAjRpVI0GTgLJ+8L2ARpiPF71hdzKa\n+l6IB6t8vZzl7W62nOXLWe57Oa6lWWs+Xpz9+j++FwIE4ODm5sb3GiI3mUw+f/7sexXw6dWrV76X\nAASAIAEAVOAeEgBABYIEAFCBIAEAVCBIAAAVCBIAQAWCBABQgSABAFQgSAAAFQgSAEAFggQAUIEg\nAQBUIEgAABUIEgBABYIEAFCBIAEAVCBIAAAVCBIAQAWCBABQgSABAFQgSAAAFQgSAEAFggQAUIEg\nAQBUIEgAABUIEgBABYIEAFCBIAEAVCBIAAAVCBIAQAWCBABQgSABAFQgSAAAFQgSAEAFggQAUIEg\nAQBUIEgAABUIEgBABYIEAFCBIAEAVCBIAAAVCBIAQAWCBABQgSABAFQgSAAAFf4faCHcJcIKEtMA\nAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "load Lshapemesh\n",
    "N = size(node,1); \n",
    "NT = size(elem,1);\n",
    "valence = accumarray(elem(:),ones(3*NT,1),[N 1]);\n",
    "valenceNew = accumarray(elem(:,1),ones(NT,1), [N 1]); % for newest vertex only\n",
    "allNodes = (1:N)';\n",
    "intGoodNode = allNodes((valence==valenceNew) & (valence==4));\n",
    "bdGoodNode = allNodes((valence==valenceNew) & (valence==2));\n",
    "figure(2); showmesh(node,elem);\n",
    "findnode(node,intGoodNode,'noindex');\n",
    "findnode(node,bdGoodNode,'noindex');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 2 -- Remove Good Nodes\n",
    "\n",
    "The star of good nodes are found using the incides matrix between\n",
    "triangles and vertices."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "t2v = sparse([1:NT,1:NT,1:NT], elem(1:NT,:), 1, NT, N);\n",
    "[ii,jj] = find(t2v(:,intGoodNode));\n",
    "nodeStar = reshape(ii,4,length(intGoodNode));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For interior good nodes, there are 3 configurations of 4 triangles in the\n",
    "star of a good node showing the following figures. Other combinations are\n",
    "ruled out by our ordering assumption: the left child is stored in a prori\n",
    "of the right child. Therefe `nodeStar(:,1)` is always in L type and\n",
    "`nodeStar(:,4)` is in R type. Here the oritentation is respect the\n",
    "direction from the new added vertex to the opposite vertex in each\n",
    "triangle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBygkOKU6MgAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNTo0MDozNh7JJbIAAB+1\nSURBVHic7d1vbFt3vcfxT9qyFQ0lXrUxKDU42i4a5HJvhsRGotE/gytaQGxVp1Ri05ZIDGnSQEWF\nB6gPmiJR8adl3aTxgDGaoj1ZtHV3F6ntVCAJV2toQVDGOlVo084wUMSk4SYU7iBZ7oNjXMexj3/H\nPn++x36/1Aet5/jr5fPN+Zw4jt2ztLQkAADStirtOwAAgEQhAQCMoJAAACZQSAAAEygkAIAJFBIA\nwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAh\nAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJ\nFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIA\nmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygk\nAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCB\nQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAA\nEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUE\nADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQ\nSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBg\nAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAA\nACZQSAAAEygkAIAJFBIAwAQKCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIAwAQK\nCQBgAoUEADCBQgIAmEAhAQBMoJAAACZQSAAAEygkAIAJFBIQL8/z0r4LiB6xxoFCssLzvLGxMba8\nw3ie19/f39/fT7KdhFhj0rO0tJT2fYA8z7tp43+WinODH/mPp3/wTKFQSPseIQJ+rLkNfbl8ry6s\nJtnOQKzx4Tuk9Pn7fceDWyVtPnDL9ntu57SrA1RiLQxvyOV7SbYzEGusKKSUVe13XhIr3hlqYhXJ\ndgRijRuFlKaV+y1WPPvqxiqSzThiTQCFlJpG+y1WPMsCYhXJZhaxJoNCSkfwfosVz6amsYpkM4hY\nE0MhpcBlv8WKZ41jrCLZTCHWJFFISXPfb7Hi2REqVpFsRhBrwiikRIXdb7HiWdBCrCJZ84g1eRRS\nclrbb7HitrUcq0jWMGJNBYWUkHb2W6y4VW3GKpI1iVjTQiElof39FituTySximSNIdYUUUixi2q/\nxYpbEmGsIlkziDVdFFK8ot1vseI2RB6rSNYAYk0dhRSjOPZbrHjaYopVJJsqYrWAQopLfPstVjw9\nscYqkk0JsRpBIcUi7v0WK56GBGIVySaOWO2gkKKXzH6LFU9WYrGKZBNErKZQSBFLcr/Fiicl4VhF\nsokgVmsopCglv99ixeOXSqwi2ZgRq0EUUmTS2m+x4nFKMVaRbGyI1SYKKRrp7rdY8XikHqtINgbE\nahaFFAEL+y1WPGpGYhXJRopYLaOQ2mVnv8WKR8dUrCLZiBCrcWvSvgPZZm2/VbXiT//gmUKhkPbd\nyaRKrLl83/TBU/6FpeJcYTg/ODKQzH0oFefOTr5QM5pk20Gs9vUsLS2lfR+yqrqN6u5Z2BscX39w\n/I+7Q31Io7ml4tz0l06z4i3wPK+/v78wlC/9/mKpOLfyCrl8b2EoPzgy4HgW4h/7Nu8ebnrNUnGu\nVLzozRbPTp5rNDq3oU8XVk9NTZFsKMSaCRRSi/z93rx7SFLAnoVaccdCctnvwlC+MJz3nvojnRSK\nH2v1JYWCRke1aZNmZjQ9renpqv80lL/j0NZcvjf4Nh2PXKXi3PTBU2cnzzmNLhQ4eLkj1qygkFrh\n73cu31vdB+2vuEshhdrvXL43t3pdl6+4u8phq/Ip3bx55XU0Pa2ZmfJfcvne0Sd3BifrcuQqFef+\ne9cJb7boPpqDlyNizRAKKbTqs61oV7xpIbWw3+K0y00l1tFRHT7s9CH79ml8vHmyTY9clVjDjibZ\npog1W3iWXTjVbTQ6qlde0d69dVpB/+qqw4c1NaXxcZWKcxN3PlH3ETZH1fvtPte/z1u2bOFpPAFa\nOGxJ2rs3gmRbOGxVRpNsMGLNHAophJo2SnLF29lv0UmBWjts+dpMtrVYq0eTbCPEmkUUkquW28jX\nzoq3v9+ikxpo57DlqyRbeSaxu+mDp1qLtXo0ya5ErBlFITnxn+HtP6Cc/IpHst+S9M5FfgWvWvVJ\nRmufW5//8OnZyXPeqaL7R52dPOc/OaX90d158GqEWLOLQmrOb6PCUN7/zsZxz77yFQ0P68yZZRe2\nsOLu+33kiIaH1durG27Ql7+sixdr53qzxcKO9XSSr+Ykw13dZP2f53mzIY5c/g40Hb2woP37NTCg\nq67SwIAOHdLiYp3RuZuvIlllJ9aKV17RrbdqeFi/+12d0d0WK4XUROW3X/1/Ou7ZmTN66CHNzi5r\nBV/YFXfc7/37NTqq2VnNz+vll3XggD796WVHLn9uqXiRlyrRilg3bXL9wEbJ+rdQ/XT85vdhtugy\n+v77tWePXnxRi4t68UV98Yv6/OfrjPZmiySboVh9i4u6+24995xmZ/X3v9cZ3W2xUkhBql+LwXHP\nDh3Sjh36yEdq16si7Iq7zH39dX3965I0Pq5SSQ8/LEk//al+8pM6c3n5rJqX2JDbeUZwsv5zdEvF\nOfcfELqMfvllPfaYJD3+uObn9a1vSdKjj+rSpTqjJXVzshmKtWL/fp1q8Ph9d8ZKITVU8zp1jns2\nMaGjR/WPfzS8QtgVd5n7y19qfl7r1mnPHvX16YEHtGGDJP3mN/XndnMnVcfqnxY4HjuCky0ULn8P\n6nJrjqN//nMtLen663XXXXrLW/S5z0nSwsKyE/nq0V2bbLZi9Z0+ra9+VZ/4RPPR3RMrhVRfTRu5\n79m3v62jR3X0aMMrhFpxx7k9Pdq2TXfdpTVrJGlhoXwG/a53NZzbPSterSZW/+FQx19AbJqs/+md\nPjjrdE/cRu/cqYUF/fa3WljQn/+sb35Tkj7wAa1f33B0FyabuVgl/fWvuvtuXX21vve98iU9PUGj\nuyRWCqkOb8VreLvv2W23aft2bd8edB33FXec+9GP6tix8iN1S0vatUt/+YtyOd12W9DcLlnxipWx\n+t7zHqcPb5qsfztNXyAq1OieHq1erVWrNDmp667T176mfF4//nGT0V2VbBZjlbRrl156Sd/9rq67\nznV0N8RKIdVqtN9yXvGmwq64+9zXXtOnPqXvfEdXXKHvf1/XXttkbjesuK9urP7Doa++Gs2IULcT\ndvS7362tW5XLqVjUzp2an28yukuSzWisR4/qscc0NqY77gg3uuNjpZCWadRGaa14qLmzs7rpJh07\npg0bdPJknZO+urfT8SuuxrH6L9Ze/XK07Si/cqDbK7s7jv7nP/XGG1pc1K236vhxnTunt75VU1M6\nfrz56I5PNruxPvKIJBWL2rFDd95ZvvCBB/TQQ81Hd3asFNJlAd8bpbXi7nOnpvSxj+kPf9C2bfr1\nr7VxY4i5nb3iAbHm8n2SovqfLn96h5yOXI6jP/tZrV2rPXvK/1y/vvxtbs2vyzQa3cHJZjpW/xWt\nf/QjHT2qp58uX3jypH7xC6fRHRwrhVQWsN9Kb8Ud5166pM98Rn/7mz75Sf3wh1q3LvTcTl3xZrH2\nFobylZdFb8fERPkGHR+JdRx9442S9MQT5afV/epXOn9ekt7/ftfRHZls1mMdH9eTT5b/TE6WL3z4\nYX3hC66jOzJWUUi+4P1WeivuOPfxx/WnP0nS8eO68kqtWVP+s7vqvSyazu28FW8aq6TC8AZJR460\nO2tmRnI+j3Yffd99uuYaeZ7e9z5t2aKbb5ak/v5lj8c2Hd1hyXZArBs3aseOy398H/+4PvShEKM7\nLFYfheS030pvxV3mPv98+S9vvqnFxct/3nwz3NxOWnHXWIfyiuLB2FA/aXAffc01evZZ3XKLLlzQ\n9LQWFrRtm06c0NVXhxvdMcl2RqxRje6YWCu6vZAc91vh92xpSUtL+q//qr087Iq7zH3kkfK4mj8P\nPhh6bmesuHusuXzf4MiA52nfPtcbX5ns2Jg8T4MjA4MjA+530nH0Bz+on/1MpZJeeEHz8zp2TO99\nbyujOyDZToq1YtWq8uiujbXa6vHyC0F3I/f9liT1/N/cG+efe62np/4747kYG9P0tAZHBuq+0eT0\nwdl6l8c+t8bavivf8eFrv3H/ga1btuVyuRZHpidUrGv7rnzHwNvPP/vSif95o7XP8NiYJiaUy/eO\nPrWz/v2ZLareqUCo0WvX6u1v1xVXhBtdZ2Jmk+3IWFsbXWdiZmOt0b2FFLKNktizuoWU/H4ryyse\nNlZJa/uuvHHrv7X2Ga58bnedvq/hXWpw5EpgdN2JWUyWWINlNNaVurSQWthvxb9nDb5DSmG/lc0V\nby1WtfoZdvzcBhy54h7daGK2kiVWF5mLta6eJf8p8d2k5f32Vd71dXxce/c6fYjjno2vPzj+x92N\n/mt8cwOUinPTXzr99A+eKTi+NFh62oxVVZ9h/3X/Nm2q8yqCnifP08yMJibkeU6fW//9GIMfLI1p\ndPDETCRLrKFkJdZGuq6Q2t9vxbZnwYUU39xgmVjxSGJV1bvF+//0/49HR7Vpk2ZmND297NklhaH8\n6FMjTW/T5cgV0+imE40nS6wtsB9rgO4qpKj2W/HsWdNCimluU8ZXPMJYVX6HjotnJ8+VinN130cx\nl+8dHBkoDOUdxzkeueIY7TLRbLLE2jLLsQbrokKKdr8Vw565FFIcc12YXfHIY63mf4b911z3f6HY\n5QBUw/3IFflox0EGkyXWNtmMtaluKaRY91sR7ZljIUU+132WtRWPO9ZItHbkSpK1ZIk1EtZidbEm\n7TuQhAT2O5fvHcyH+wW6zM2t/AqekRXPxGErE0wlS6xRMRWro85/pYYM7Xch7TvQlJ1fCyfWaBlJ\nllijZSRWdx3+kF1W9ntw8tztu05IKuV7z44MzBh+HEAGHgrISqy54tztu04UZouSvKH8M4e2lsK8\n8Wjy0k2WWGOS+hesu04upKzsd+FU8d47J6svOfLkiGf7Pqe44lmJVdK9OyYLVc868YbyR6J40mOs\n0kqWWGOVlU7q5ELq6ekpDOX9V8u2bO+p329e/mS5iaH8EfN3u1ScK525NDU1lfCKZyVWSVMHZ2su\n2bJ7KJV7EkoqyRJr3NL6gg2lk3+GVCgUHN9WK13eiksKid+HDMlKrAXnCyFihaTOfpbd1NTU9ntu\nH9z978YX/eLInG55tPqSV3cPbbb9wEWpOPffO08ePnw4+bOtrMQqyTv1+5rHdgq7hwvp3R8XaSVL\nrLFK8Qs2lA7/DunpHzwz/aXTpeJc2vclSCnfO/2vb/lL+d5nDm21/wMkf7k3t/x+GG3ISqySnjm0\n9ey/npGflZ80pJUsscYn3S/YUDr5Z0g+z/O233P75gO32D/zmlh/cDTkL8Ymz8hyZyhW7+Apz/Zv\nUPosJEuskbMQq7tO/g7Jl6EzLy/tO9CUneUm1mgZSZZYo2UkVned/DOkCn/LYz3zKhXnzk6+UPl7\nYTifzKsnJDnX2nInEGuXMJUssUbFVKyOOv8hu4rIHw3wX+TUmy36r3O68gq5fG9hKD84MhDHi6tG\nONeF2eWO9UGeSPq+5VfhTOZUw2ayxNomm7E21UWFpEi3vFScmz546uzkucolhULQ20DccWhr06GO\nbz8R+dymjC93tAevyPs+7PsUJHmqYTlZYm2Z5ViDdVchKaItr7wpUaUMVubueZqe1sxM+S+5fO/o\nkzuDh7q8QV8cc4NlYrmjOnjF0ffu7+SW8KmG/WSJtQX2Yw3QdYWktre80gqjozp82OlD9u3T+Hjz\nbmj6FuYxzQ2QoeVu/+AVU987vtd1wqcaWUmWWEPJSqyNdGMhqY0tb6EVfC7dEFBIsc5tJHPL3c7B\nK76+b3rkSv5UI1vJEqujbMVa1+rx8fG070MKcrnc1i3bvnH/gXd8+Nq1fVc6flTLrSBp82b19OjE\n/7xx/tmXbtz6b3WHTh+crbvfcc+tK4vL3VqsavUz7Pi59d/Vt9GPB2Id3XBippIlVheZi7WuLi0k\ntbTlJ/ZOnX/2pRZawVfZtv+be+PGrTesvEKjQop77krZXe4WYo277wOOXMmfamQ0WWINltFYV+r8\nX4wNEOq38M5OnvN/MtlaK/j27tXmzTo7ec5/03EXyc/N+nKH/eXK6YOnWjt2+Pbu1fh4+WfXYT82\n4dGZTpZYG8l0rDW6upAUZsv9Q/noqOstv/KKbr1Vw8P63e+WXe7vjDfrWkiOcxcWtH+/BgZ01VUa\nGNChQ1pcbGVuZyy3e6wt9P3KZFs4zwg1+sgRDQ+rt1c33KAvf1kXL17+T+6jOyDZOGIN+MIh1uR1\neyHJecv9Q/mmTU63ubiou+/Wc89pdlZ///uy/+TfQvXTQIM5zr3/fu3Zoxdf1OKiXnxRX/yiPv/5\n0HM7abldYw15ntEo2bDnGe6j9+/X6KhmZzU/r5df1oED+vSnl51tuIzumGQjjzX4C4dYE0YhSW5b\n7v8nxyPX/v061eC7bf/V30vFOcdHHlzmvvyyHntMkh5/XPPz+ta3JOnRR3XpUoi5nbfcLrGGOs9Q\n42TDnmc4jn79dX3965I0Pq5SSQ8/LEk//al+8pMQozss2QhjbfqFQ6wJo5DKgrfcXwvHNjp9Wl/9\nqj7xiUaDyuc+peLF+tcIP/fnP9fSkq6/Xnfdpbe8RZ/7nCQtLCx7ECB4bkcutxwOXqHOMwKSDXue\n4Tj6l7/U/LzWrdOePerr0wMPaMMGSfrNb1xHd2SyUcXa9AuHWBNGIV0WsOX+t+Eub23117/q7rt1\n9dX63vfKl/T01F7HX6HpFW+EvJLj3J07tbCg3/5WCwv685/1zW9K0gc+oPXrneZ26nL7AmINdZ4R\nnGyo8wz30T092rZNd92lNWskaWGhfPL+rnc5je7gZCOJtekXDrEmjEJaJvjM6z3vaX4Lu3bppZf0\n3e/quusaXse/HfdffGs6t6dHq1dr1SpNTuq66/S1rymf149/7DS3g5e7olGs7ucZckjW/TzDffRH\nP6pjx8oP6Swtadcu/eUvyuV0223NR3d8su3H6vKFQ6xJopBq1d1y/++vvtrkY48e1WOPaWxMd9wR\ndLWmtxN2bsW7362tW5XLqVjUzp2an28yt7OXu1rAqYbLeYZLsmHPMxxH+157TZ/6lL7zHV1xhb7/\nfV17bZPRXZJsm7FWBHzhEGuSKKQ6Vm755t1D0rKXQazrkUckqVjUjh26887yhQ88oIceWnY1/3Zc\nXtnXce4//6k33tDiom69VceP69w5vfWtmprS8eNBczt+uWusjNW9712SdT9vCDVa0uysbrpJx45p\nwwadPKnt22uvUHM7XZVsO7G6fOEQa5IopPpqtjyX75PkeU0+yn9dwB/9SEeP6umnyxeePKlf/GLZ\n1crFMNS8kBznfvazWrtWe/aU/7l+ffnc6syZhnO7YblXqonVse/llqz7eUao0VNT+tjH9Ic/aNs2\n/frX2rixznWqR3dhsi3H6vKFQ6xJopAaqt5y/w1L/JfjDTA+riefLP+ZnCxf+PDD+sIXLl9nYkKS\ncvlelwcBHOfeeKMkPfFE+dlBv/qVzp+XpPe/v/7cLlnuupbH6tT3ckvW/TxDzqcaly7pM5/R3/6m\nT35SP/yh1q2rf7XK6K5NtrVYXb5wiDVJFFKQ6i0vDG+QdORI0PU3btSOHZf/+D7+cX3oQ5evMzMj\nOe+3JJe5992na66R5+l979OWLbr5Zknq71/2IEBlbvcsdyOVWCW59L0ckg11niHnU43HH9ef/iRJ\nx4/ryiu1Zk35z+6qF4WvjJbUzcm2EGvTLxxiTRiF1ERly/1zH5dvxoOFegRA/6qu4LnXXKNnn9Ut\nt+jCBU1Pa2FB27bpxAldfXXt3Fy+t3uWO0BVrL1q1vcuwp5nyO1U4/nny395800tLl7+8+abtaNz\nG/pINmysTb9wiDVhFFJz/pafPXTeP/fZt8/po1at0tKSlpb03vdevnBsTJ6nwZGBwZEBx+m5fN/g\nyEDTuR/8oH72M5VKeuEFzc/r2LE6cwtDee+pC92z3MH8WEtnLinkeUbdZMOeZ8jtVOORR8qzav48\n+GDtaF1YTbIKH2vwFw6xJoxCclIoFKampnRhtaTxcddOqjE2pokJ5fK9dxza6v5RuXzv5t3DuXyv\ny9y+Pg0M6G1vqzNX6rrlbsqPtVAouJ9n1NXCeYacTzVcRhcKBZKtaCHWRl84xJowCslVZcvVUidV\n2mjX6fvCjq68g2TLc6VuXG4XlVgTPs9QyFONgNEkuxKxZhSFFELLndROG/la6yTayEU7B682k23/\nVINkGyHWLKKQwqnppP7+qgfElvOfabNvn/r7220jX/WWu88VbeSgtYNX++cZau9Ug2SDEWvm9Cz5\nv/KHMDzP27JlS2nx9covh/uvXjU6qk2bNDOj6ellP9UsDOVHnxpperPj6w+O/3F38HUq747sMjeX\n782tXte1yx2WH6vnef7rWm7aVOdVMj1PnqeZGU1MyPOcDlv++37WfXP6ilJxbuLOJ0rFOffR3XzY\nCoVYM4RCapG/5YO7bvROFUvFubrvo5XL9w6ODBSG8o7P0nEpJJVfkf7i2clzwXNz+b6zh85383K3\noHLw8v8ZyXmGy5FLIU81uvywFRaxZgWF1DrP87bfc/vmA7f4r33gzRb91/r1fwei6aau5FhI1RrN\n5bdfW+YfvPTORTV4v86w5xmORy65nWoUhvKlM5dINixizQQKqS3VndT+rbVQSHXRRm2qxCqp/fMM\n9yNXtbqnGiTbDmK1j0JqV4SdFEkhsdyRiDDW1o5cK5Fs+4jVOJ5l167g9/RLGMsdFVOximQjQqzG\nUUgRMLLlLHe0jMQqko0UsVpGIUUj9S1nueOQeqwi2RgQq1kUUmRS3HKWOz7pHrxINibEahOFFKVU\ntpzljltaBy+SjRWxGkQhRSzhLWe5k5H8wYtkE0Cs1lBI0Utsy1nuJCV58CLZxBCrKRRSLBLYcpY7\neckcvEg2YcRqB4UUl1i3nOVOS9wHL5JNBbEaQSHFKKYtZ7nTFd/Bi2RTRKwWUEjxinzLWW4L4jh4\nkWzqiDV1FFLsItxyltuOaA9eJGsEsaaLQkpCJFvOclsT1cGLZE0h1hRRSAlpc8tZbpvaP3iRrEHE\nmhYKKTktbznLbVk7By+SNYtYU0EhJaqFLWe57Wvt4EWyxhFr8iikpIXacpY7K8IevEg2E4g1YRRS\nChy3nOXOFveDF8lmCLEmiUJKR9MtZ7mzyOXgRbKZQ6yJoZBSE7DlLHd2BR+8SDajiDUZFFKa6m45\ny511jQ5eJJtpxJoACillNVvOcneGlQcvku0AxBq3nqWlpbTvA+R53vZ7bj/7v88XCgWWu2P4seqd\ni5JKZy6RbGcg1vhQSFZ4nrdv3757772X5e4knudt2bJFEoetTkKsMaGQgHh5niepUCikfD8QKWKN\nA4UEADCBJzUAAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBI\nAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAC\nhQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAA\nJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJ\nAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEyg\nkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDA\nBAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEB\nAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkU\nEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCY\nQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQA\ngAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFC\nAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAAT\nKCQAgAkUEgDABAoJAGAChQQAMIFCAgCYQCEBAEygkAAAJlBIAAATKCQAgAkUEgDABAoJAGAChQQA\nMIFCAgCYQCEBAEygkAAAJlBIAAAT/h8JtkOhtZENsQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "figure(1);\n",
    "node1 = [1,0; 0,1; -1,0; 0,-1; 0,0];\n",
    "elem1 = [5 3 2; 5 1 2; 5 4 1; 5 3 4];\n",
    "subplot(1,3,1); showmesh(node1,elem1); findelem(node1,elem1);\n",
    "findnode(node1,5,'noindex');\n",
    "elem1 = [5 3 2; 5 4 1; 5 1 2; 5 3 4];\n",
    "subplot(1,3,2); showmesh(node1,elem1); findelem(node1,elem1);\n",
    "findnode(node1,5,'noindex');\n",
    "elem1 = [5 3 2; 5 4 1; 5 3 4; 5 1 2];\n",
    "subplot(1,3,3); showmesh(node1,elem1); findelem(node1,elem1);\n",
    "findnode(node1,5,'noindex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBysEKOZJOQAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNTo0MzowNMlfuLMAACAA\nSURBVHic7d1/dJTVnfjxTyDyY7VhRGtbYHQQ2orpno2tlSYWMrC1RawCgkEpYiJi3R53D57gHx57\nCrFbzp5dEPBb9/gFK8Ha3UNQfthz8AfdMqmVCEpF11iqAo878q1ClSEhBTFhvn/c7BAmk5nJZJ7n\n3vs879fhj8kzl+HmerzvPDNPZkqSyaQAAKDbIN0TAABAhCABAAxBkAAARiBIAAAjECQAgBEIEgDA\nCAQJAGAEggQAMAJBAgAYgSABAIxAkAAARiBIAAAjECQAgBEIEgDACAQJAGAEggQAMAJBAgAYgSAB\nAIxAkAAARiBIAAAjECRAp/Xr199444033njjjh07dM8F0KxU9wSA4Nq7d+8999xz+vRpEZk5c6bu\n6QCacYYE6NHW1jZ37lxVIwBCkABdFi1adODAAd2zAAxCkAANHnvssaampmHDhi1YsED3XABT8BoS\n4LU33njjvvvuE5E1a9acOHFC93QAU3CGBLjLcZyeX544caKmpubUqVO33nrr3XffrWlSgIkIEuAW\nx3GmTJkyduzYkpKSKVOmqDLdc88977zzzvjx49euXat7goBZeMoOcEtdXV0sFlO3Y7FYXV3d7bff\n/qtf/Wro0KFNTU2f+9zntM4OME5JMpnUPQfAn0pKStKO3HLLLZs2bRo+fPjo0aPVkWPHjn388cci\ncskll5SVla1Zs2b69OleTxQwA2dIgCvSXjpSOjo6ROTkyZPvvfde2l1Hjhw5cuRIe3u7B3MDzESQ\nAFdEIpFoNJp6yk5EotHozJkzx4wZ03PYm2+++corr4hIdXX1V7/61XHjxnk8T8AcBAlwheM4+w78\noaKmfF9Tq4hEKsOJrk+uu+66RYsW9Rz28MMPqyDdfvvtCxcu1DNXwAwECSg+x3Gumvx3M1dNi1SF\nQ+EyEYnWVyXibbMWzNjy5LZIJKJ7goCJuOwbKLKeNep5PBQui66YOGvBjIwvLwHgKjugmHrXKLZy\nl4hE66vUl4l4W2zJbs6TgN44QwKKpq9zo544TwL6QpCA4sinRgpNAjIiSEAR5F8jhSYBvREkYKD6\nWyOFJgFpCBIwIIXVSKFJQE8ECSjcQGqk0CQghSABBRp4jRSaBCgECShEsWqk0CRACBJQgOLWSKFJ\nAEEC+seNGik0CQFHkIB+cK9GCk1CkBEkIF9u10ihSQgsggTkxZsaKTQJwUSQgNy8rJFCkxBABAnI\nwfsaKTQJQUOQgGx01UihSQgUggT0SW+NFJqE4CBIQGYm1EihSQgIggRkYE6NFJqEICBIQDrTaqTQ\nJPgeQQLOYWaNFJoEfyNIwFkm10ihSfAxggR0M79GCk2CXxEkQMSeGik0Cb5EkADLaqTQJPgPQULQ\n2VgjhSbBZwgSAs3eGik0CX5CkBBcttdIoUnwDYKEgPJHjRSaBH8gSAgiP9VIoUnwAYKEwPFfjRSa\nBNsRJASLX2uk0CRYjSAhQPxdI4UmwV4ECUERhBopNAmWIkgIhODUSKFJsBFBgv8FrUYKTYJ1CBJ8\nLpg1UmgS7EKQ4GdBrpFCk2ARggTfokYKTYItCBL8iRr1RJNgBYIEH6JGvdEkmI8gwW+oUV9oEgxH\nkOAr1Cg7mgSTEST4BzXKB02CsQgSfIIa5Y8mwUwECX5AjfqLJsFABAnWo0aFoUkwDUGC3ajRQNAk\nGIUgwWLUaOBoEsxBkGAralQsNAmGIEiwEjUqLpoEExAk2IcauYEmQTuCBMtQI/fQJOhFkGATauQ2\nmgSNCBKsQY28QZOgC0GCHaiRl2gStCBIsAA18h5NgvcIEkxHjXShSfAYQYLRqJFeNAleIkgwFzUy\nAU2CZwgSDEWNzEGT4A2CBBNRI9PQJHiAIME41MhMNAluI0gwCzUyGU2CqwgSDEKNzEeT4B6CBFNQ\nI1vQJLiEIMEI1MguNAluIEjQjxrZiCah6AgSNKNG9qJJKC6CBJ2oke1oEoqIIEEbauQPNAnFQpCg\nBzXyE5qEoiBI0IAa+Q9NwsARJHiNGvkVTcIAESR4ihr5G03CQBAkeIcaBQFNQsEIEjxCjYKDJqEw\nBAleoEZBQ5NQAIIE11GjYKJJ6C+CBHdRoyCjSegXggQXUSPQJOSPIMEt1AgKTUKeCBJcQY3QE01C\nPggSio8aoTeahJwIEoqMGqEvNAnZESQUEzVCdjQJWRAkFA01Qj5oEvpCkFAc1Aj5o0nIiCChCKgR\n+osmoTeChIGiRigMTUIagoQBoUYYCJqEnggSCkeNMHA0CSkECQWiRigWmgSFIKEQ1AjFRZMgIqW6\nJwD7aKxRIt62r+mt1O1IVbiiptzjOfiSCQubatKWJ7dFIhGP/3WYoCSZTOqeA2zifY0S8bZE/LjT\nEt/X1JqIt/UeEAqXRSrDFTXlZp6uxVbuEpFofZXuiaQzc2ET8bbYkt00KZgIEvpBS41iK3fta2pN\nHYlEpLZWqquluVliMYnFzg6OVIZnrp4WCpd5M7c8mRkkkxeWJgUWQUK+tNRo6+LnnZZ4aq+MRnvP\nSmIxaW7uvhEKl9U+PdeoJhkYJPMXliYFE0FCXjTWqLZW1q/P6680NMiyZcY1ybQg2bKwNCmAuMoO\nuVlRIxFZulSWLZNEvK1xzsaML4rAooXlursAIkjIwZYaKTQpC+sWliYFDUFCNlqu8I6t3FXYpqmk\ntk71XBlSbFxYmhQoBAl90lKjfU2t6tKvnJvmAw9IVZXs2ZPhrqVLJRqVfU2tzq64C3O0Uj4Lu2GD\nVFVJWZmMHy/33y/Hj6cP0LKwNCk4CBIy0/Xbr2qnq63NMWzPHlmzRlpaMmyairpszGkhSN1yLuzy\n5VJbKy0t0t4uBw7IihVy003S1ZU+TMvC0qSAIEjIQON7Maidrrq6zwGrV8vs2TJpkpw8me1x1CP0\n/D2bgMu+sJ98Iv/yLyIiy5ZJIiGPPCIi8rvfyW9/mz5S18LSpCAgSEin933q1AvmWX6Qb2yUzZvl\n9Okcj6MuFU7E27i0Qcm+sH/4g7S3y8iR8uCDMmKE3HuvjBkjIvLf/50+UuPC0iTfI0g4h94aqZ+7\nsz9f9/DDsnmzbN6c46Eike4nlxLxPp7UC5KcC1tSItdfLz/4gZSWioh0dkpHh4jI6NHpI/UuLE3y\nN4KEs7S/h7d6nSP770FOnSqzZsmsWbkfTe2bsZUtRZiZ5XIu7N//vWzf3v1MXTIpixfLsWMSCsnU\nqRkG611YmuRjBAndtNco5bLLivk45rxlg3b5LOzRo/L978u//7sMGSJPPCGf/3yfj6NxYWmSXxEk\niBhTI/WyxPvvF+fRivU4PpDnwra0yFVXyfbtMmaM7NjR52moCQtLk3yJIMGUGolItL5S5Jz3mR4I\n9TjavykT5LOwO3fKd74jhw/L9dfLG2/I5Ml9jjRkYWmS/xCkoDOnRiISCo8QkWJtL937ZqX+70u7\nnAvb0SHz5slf/yo33CC//rWMHJnt0cxZWJrkMwQp0IyqkfzvJ8KpzzsYoMbG7gfkNSTJY2Gfeko+\n/FBE5LnnZOhQKS3t/lNfnz7StIWlSX5CkILLtBopkaoxIrJhw0Afp7lZxIyf4g2RfWHffLP7xpkz\n0tV19s+ZM+kjDVxYmuQbBCmgzKyR/O9Ol88ZUjIpyaRcd13mew15ncMc2Rf20Ue71zPtz6pV6SPN\nXFia5A8EKYiMrZGIhMIjKmrKHUcaGgp/kLo6cRypqCmvqCkv3tTs5vuFpUk+QJACx+QaidpW6qtC\n4bJlywrcOuvqpLFRQuGymaunFXt2FgvCwtIk2xGkYDG8Rkrqo7IL2DpTm+bi3YvcmZ3FgrCwNMlq\nBClArKiRUtjWacumqVEQFpYm2askmUzqngO8YFGNUlIfmK3e0LO6OsPbgzqOOI40N0tjoziOiZum\n+nzVaH2V7omc5Y+FzS4Rb4st2b3lyW2R7O+NCJMQpECwsUZKIt62dfHzqY+DU3tLba1UV0tzs8Ri\n51w2FqkM1z5T4/0kszMwSOKLhc2JJlmHIPmfvTVSEvG2RPz4vqbWRLwt4weVhsJlFTXlkcqwmd+g\nmUES+xc2HzTJLgTJ52yvURq1daoPU1DvFGDgRp/G2CD1ZOPC5okmWYQg+ZnPamQpK4LkbzTJFlxl\n51vUCFC47s4WBMmfqBHQE02yAkHyIWoE9EaTzEeQ/IYaAX2hSYYjSL5CjYDsaJLJCJJ/UCMgHzTJ\nWATJJ6gRkD+aZCaC5AfUCOgvmmQggmQ9agQUhiaZhiDZjRoBA0GTjEKQLEaNgIGjSeYgSLaiRkCx\n0CRDECQrUSOguGiSCQiSfagR4AaapB1Bsgw1AtxDk/QiSDahRoDbaJJGBMka1AjwBk3ShSDZgRoB\nXqJJWhAkC1AjwHs0yXsEyXTUCNCFJnmMIBmNGgF60SQvESRzUSPABDTJMwTJUNQIMAdN8gZBMhE1\nAkxDkzxAkIxDjQAz0SS3ESSzUCPAZDTJVQTJINQIMB9Ncg9BMgU1AmxBk1xCkIxAjQC70CQ3ECT9\nqBFgI5pUdARJM2oE2IsmFRdB0okaAbajSUVEkDR45plnbr755jFjxlw+7vIh5w85+s7HZzrP6J4U\nNPiw9ch/3LHlP+7YcvTdT3TPBYXrb5O6urpWrVo1efLkkSNHjhkzZvr06S+++KLLc7RDSTKZ1D2H\nAEkmkz/60Y8ee+yxtOOjK75Yt2Vu6dBSLbOCq2Ird4lItL6q913b7nvh9Y1vicid22699JujvZ4Z\niioRb4st2b3lyW2RSCTLsFOnTk2dOrWlpSXt+OLFi1etWuXi/GzAGZKnNm3a1F2jEglfPSp89Sh1\n/PC+D198qFnnzOCtjo9PvvhQs6oR/CHP86Sf/vSnqkbDhg2rq6u76aabhgwZIiKrV6/eunWrN1M1\nFkHy1NKlS9WNKUuqFj5728Jnb7v559PVkdZn/6RvXvDOm5v/+H+ufeLhrz+267HXdM8FRZazSclk\n8uc//7mIlJaWHjhw4Iknnti2bdtPfvITde/27ds9m6qZCJJ3Ojo69u/fLyJlX/pc9X2V6mD5jV8Z\nVjZURDo+Ptn+UYfO+cETHUf/+vGhY12f8aqhP2Vv0rvvvtvW1iYikyZNGjWq+wmSGTNmqBunTp3y\napqG4kUL7+zevbukpKR02ODLvjUmdTB5Jtn5aZeIlA4tPf+i4fpmB49cOnH0dx6cpG7vferNY+8f\n1zsfFF2qSb1fTxo7duxHH30kIsOHd//PfubMmbVr16rb1113nbczNQ5BcldDQ4PjONXV1dFodHbt\nrDs23ZJ2hfdLj+zu/LRTRMLfHDWolBNWv0nE25xdH4iIUxlX/+lHV3xxdMUX1b3v/fYQQfKlnk3a\nsGGD2gRqa2vPO++8Sy65RI354x//+MMf/nD//v1Hjx4dPHjwj3/843nz5umdtnYEyUVTpkyJxWIi\n0tjYKCK1T9f0rFHXZ12/+dlLLWv3ikjpkMHT/3mqnlnCNc6ueOOcJnW7cU585uppFTXleqcEz6gm\njR07Vn3Z2NjY0NBw6NCh1ID29vaXXnpJ3b7gggsuu+yyQYOC/iNp0L9/9zQ2NqoapTgt8dTtD98+\n+n+/95Sq0fDQsPn/OefzX7nI4xnCbakaKer6bwRHz//lRcRxnIaGhtSX48eP/+Uvf9nQ0BAOh48f\nP37nnXcuWrTI8zkaJgl31NbW6v5vC8As8+fPTyQSiUSis7MztVd88MEHJSUlInL++eefPHlS466l\nHU/ZuaW6ulo9U5cyc/W0ilvKn18We2XdXhEpHVr63Z9MvqbuKj3zg/uWjVqZfuT/1ff8snH2Rqfl\nA+EXY31qX1Pr1sXP9zzyySefhEIhEdm2bdtNN92kDo4ePXrUqFGHDx/u6OjYu3fvtddeq2GuZuAp\nO7dEo9G0C2wileFXHt+rajTkb867Y9Mt1MjfZq6e1vPL2qdrdM0EWoTGlPX8MhKJzJ07V91+5pln\nUsdfffXVw4cPq9sTJkzwbHoG4gzJLZFIZOfOnbFYrLm5ORKJ3HHHHTfe9v1333pX3XvJhItbf/1O\n66/f6flXJi686sJLR+iYLFxRUVMeqQyrn5Gj9ZW8hW6gJOJtsSV71q9fLyJqE1i6dOnhw4fLysra\n2tqefPLJkSNHzp49+4033nj00UfVX7nmmmtGjhypddaa8V523lm3bt3dd9+dZcCdW2+99Bqet/Gb\nLO9lx1N2fpWIt22du2P9+vXRaDTtrscffzzjxQsXXnjha6+9dvnll3sxP1PxlJ13jh49qnsKAFyX\npUYictddd+3YsaO8/OwvAAwePHj+/Pl79+4NeI2EMyTvOY4za8GM6IqJoXBZ7tGwX5YzJPhP9hqd\nMzKR+NOf/nTBBReMHz9+6NChnszOdJwheS0SiWx5cltsye5EvE33XAAUU/41EpFQKDRx4sTy8nJq\nlEKQNKBJgP/0q0bIiCDpQZMAP6FGRUGQtKFJgD9Qo2IhSDrRJMB21KiICJJmNAmwFzUqLoKkH00C\nbESNio4gGYEmAXahRm4gSKagSYAtqJFLCJJBaBJgPmrkHoJkFpoEmIwauYogGYcmAWaiRm4jSCai\nSYBpqJEHCJKhaBJgDmrkDYJkLpoEmIAaeYYgGY0mAXpRIy8RJNPRJEAXauQxgmQBmgR4jxp5jyDZ\ngSYBXqJGWhAka9AkwBvUSBeCZBOaBLiNGmlEkCxDkwD3UCO9CJJ9aBLgBmqkHUGyEk0CiosamYAg\n2YomAcVCjQxBkCxGk4CBo0bmIEh2o0nAQFAjoxAk69EkoDDUyDQEyQ9oEtBf1MhABMknaBKQP2pk\nJoLkHzQJyAc1MhZB8hWaBGRHjUxGkPyGJgF9oUaGI0g+RJOA3qiR+QiSP9EkoCdqZAWC5Fs0CVCo\nkS1KdU8ALlJNmrVgRnTFxFC4TPd0iiARb9vX9FbqdqQqXFFTrndK/uDjhaVGFilJJpO65wB3OY5j\ndZMS8bZE/LjTEt/X1JrxbC8ULotUhitqyiNVYe+nl1Ns5S4RidZX6Z5IOtsXNh/UyC4EKRDsbVIi\n3hZbuWtfU2vqSCQitbVSXS3NzRKLSSx2dnCkMjxz9TTTvkczg+SDhc2JGlmHIAWFjU1KxNu2Ln7e\naYmn9sreG4vjSCwmzc3dN0Lhstqn5xr1PRoYJH8sbHbUyEYEKUDsalJq06ytlfXr8/orDQ2ybJlx\nW6dpQfLNwmZBjSzFVXYBYtF1dwVsmiKydKksWyaJeFvjnI3mf49aBGFhqZG9CFKwWNGkwjZNxa6t\n02NBWFhqZDWCFDjmNym2cldhm6aS2jrVc2VI8f3CUiPbEaQgMrlJ+5pa1aVfhW2aytKlEo3KvqZW\nZ1e8aDOznO8Xlhr5AEEKKGObpHa62to+B3R2yvLlUl4u558v5eWyerV0dWUYpjYlp8W4fVOXnAub\ncuiQfPvbUlUl//M/Ge41c2GpkT8QpOAys0lqp6uu7nPAP/yDPPigvP22dHXJ22/LfffJP/5jhmHq\nEXr+nk3A5VxYpatL5s+Xl1+WlhY5eTLDAAMXlhr5BkEKNAObpGbS1w/yBw7IL34hIvLUU9LeLv/2\nbyIi69ZJR0f6yEik+9HM+db0yr6wKcuXy66srxCZtrDUyE8IUtAZ1ST1c3eWTfPVVyWZlHHj5Ac/\nkPPOk7vvFhHp7JTjx9NHRiLdTy4l4r3uC56cC6vs3i0PPSTTp2cbY9TCUiOfIUgwqEnqdQ71M3hG\nc+dKZ6e88450dsqRI/Kv/yoi8rd/K6NGZRis9qjYyhYXZmqZnAsrIidOyPz5cuGF8vjj3UdKSjKP\nNGRhqZH/ECSImNQkEbnssj7vKimRwYNl0CBpapIvfEF+9jMJh+W//ivb41jxzgLeyLKwIrJ4sbz3\nnqxdK1/4Ql6Po3dhqZEvESR0M6FJ6p9+//3cIy+9VKZNk1BI4nGZO1fa2zOMyedxAiLnwm7eLL/4\nhdTVycyZuR9N+8JSI78iSDhLe5Oi9ZUi57zPdJrPPpNPP5WuLvn2t+W556S1VYYPl5075bnnMgxW\nj2PvRycUUc6FffRREZF4XGbPljlzug/ee6+sWZNhsN6FpUY+RpBwDr1NCoVHiIjj9Dngrrtk2DB5\n8MHuL0eN6n76aM+eDIO7981KgpR7YdV7LP/mN7J5s2zZ0n1wxw557bUMgzUuLDXyN4KEdBqbpD4R\nTn3eQUZXXCEisnFj92V1r78u+/eLiFx5ZfrIxsbuB+Q1JMljYZctk6ef7v7T1NR98JFH5J/+KX2k\nxoWlRr5HkJCBxiZFqsaIyIYNme9dtEguvlgcRyZMkClT5JprRETGjpVZs9JHNjeLcHrUQ/aFnTxZ\nZs8++0f53vfkm99MH6lrYalREBAkZKarSWqn6+sH+YsvlhdekIkT5c9/llhMOjvl+uvl+eflwgvT\nR/ICUprsC5s/LQtLjQKCIKFPWpoUCo+oqCl3HGloyDzg61+XV16RRELeekva22X7dvnKV9LH1NWJ\n40hFTXlFTbnbE7ZFzoVNGTRIkklJJk1ZWGoUHAQJ2XjfpFC4LFpfFQqXLVuWbescMULKy+WCCzLc\nVVcnjY0SCpfNXD3NvXlaJ8+FzULLwlKjQCFIyEFLk9RHZRewdaY2zcW7F7kzO4tZt7DUKGgIEnKz\npUnUKCeLFpYaBVBJUv0CApCL4zizFsyIrpjo2fW+qQ/MVm/oWV2d4e1BHUccR5qbpbFRHMfEGqnP\nV43WV+meyFnmLyw1CiaChH7Q0qSti59PfRycenvQ2lqprpbmZonFzrlsLFIZrn2mxpuJ5c/AIInZ\nC0uNAosgoX+0NCkRP76vqTURb8v4QaWhcFlFTXmkMmzmRd5mBklMXVhqFGQECf3mfZNS1NapPkxB\n/esGbvRpjA1ST4YsLDUKOIKEQmhsknWsCJIJqBG4yg6F0P6+4PAZagQhSCgYTUKxUCMoBAmFo0kY\nOGqEFIKEAaFJGAhqhJ4IEgaKJqEw1AhpCBKKgCahv6gReiNIKA6ahPxRI2REkFA0NAn5oEboC0FC\nMdEkZEeNkAVBQpHRJPSFGiE7goTio0nojRohJ4IEV9Ak9ESNkA+CBLfQJCjUCHkiSHARTQI1Qv4I\nEtxFk4KMGqFfCBJcR5OCiRqhvwgSvECTgoYaoQAECR6hScFBjVAYggTv0KQgoEYoGEGCp2iSv1Ej\nDARBgtdokl9RIwwQQYIGNMl/qBEGjiBBD5rkJ9QIRUGQoA1N8gdqhGIhSNCJJtmOGqGICBI0o0n2\nokYoLoIE/WiSjagRio4gwQg0yS7UCG4gSDAFTbIFNYJLCBIMQpPMR43gHoIEs9Akk1EjuIogwTg0\nyUzUCG4jSDARTTINNYIHCBIMRZPMQY3gDYIEc9EkE1AjeIYgwWg0SS9qBC8RJJiOJulCjeAxggQL\n0CTvUSN4jyDBDjTJS9QIWhAkWIMmeYMaQReCBJvQJLdRI2hEkGAZmuQeagS9CBLsQ5PcQI2gHUGC\nlWhScVEjmIAgwVY0qVioEQxBkGAxmjRw1AjmIEiwG00aCGoEoxAkWI8mFYYawTQECX5Ak/qLGsFA\nBAk+QZPyR41gJoIE/6BJ+aBGMBZBgq/QpOyoEUxGkOA3NKkv1AiGI0jwIZrUGzWC+QgS/Ikm9USN\nYAWCBN+iSQo1gi0IEvyMJlEjWIQgweeC3CRqBLsQJPhfMJtEjWAdgoRACFqTqBFsRJAQFMFpEjWC\npQgSAiQITaJGsBdBQrD4u0nUCFYjSAgcvzaJGsF2BAlB5L8mUSP4AEFCQPmpSdQI/kCQEFz+aBI1\ngm8QJASa7U2iRvATgoSgs7dJ1Ag+Q5AAK5tEjeA/BAkQsa1J1Ai+RJCAbrY0iRrBrwgScJb5TaJG\n8DGCBJzD5CZRI/gbQQLSmdkkagTfI0hABqY1iRohCAgSkJk5TaJGCAiCBPTJhCZRIwQHQQKy0dsk\naoRAIUhADrqaRI0QNAQJyM37JlEjBBBBAvLiZZOoEYKJIAH58qZJ1AiBRZCAfnC7SdQIQUaQgP5x\nr0nUCAFHkIB+c6NJ1AggSEAhitskagQIQQIKVqwmUSNAIUhA4QbeJGoEpJQkk0ndcwDs5jjOrAUz\noismhsJlve+NrdwlItH6qtSRHf/8u6Pvfiwin53s/HDv0YqKiosuukhE5s2bd9ttt3k1a8A4pbon\nAFhPnSdlaVKa1//zrb8eO5n6sqWlRd34xje+4dYUARsQJKAI8m9S+0cnVI1KS0uvvPLKUCjU80Hc\nnidgMoIEFEeeTfro7b+oGw899NADDzzg1ewAC3BRA1A0+Vzj4LTE1Y158+Z5NS/ADgQJKKbsTUrE\n215f3yoiZWVl+/fvr6qqGjFixBVXXHHXXXcdOXLE88kCZuEqO6D4el53l7rKTl3hPWjQoIMHDw4e\nPLirq6vnX/nSl7700ksvjRs3TtOUAf0IEuAK1aRE5zH1HF2kMix/Hrxu3bobbrjh9OnTIrJw4cLp\n06f/5S9/Wb58+fvvvy8is2bN2rx5s+Z5A/oQJMAtU6ZMicViqS+j0eizzz7b2NgoIqNHj7755pvV\n8UOHDk2YMOHTTz8VkWPHjvW87g4IFIIEuKWkpCTtSF//u02aNOn3v/+9iLS0tHzrW99yfWaAkbjs\nG3CF4zi9D27btu348eMicvXVV1955ZWp46lQDR061JPZASbiDAlwS++n7ObMmXPvvfeKyLXXXqtO\niUTk4MGDEyZMOH369LBhw9rb20tL+TERAcVl34Bb1q9fX1tbq25Ho9GdO3dOnTp10KBBIvLyyy8v\nWLDghRde2LRp03e/+111mcOdd95JjRBknCEB7nIcp+d7Aq1YseL+++/vPexrX/vanj17hg8f7t3M\nAMNwhgS4K+0d6pYsWbJx48aeLyANGzasvr5+9+7d1AgBxxkSoMeRI0cOHjwYCoW+/OUvDx48WPd0\nAP0IEgDACDxlBwAwAkECABiBIAEAjECQAABGIEgAACMQJACAEQgSAMAIBAkA3HMoEwAAAJ9JREFU\nYASCBAAwAkECABiBIAEAjECQAABGIEgAACMQJACAEQgSAMAIBAkAYASCBAAwAkECABiBIAEAjECQ\nAABGIEgAACMQJACAEQgSAMAIBAkAYASCBAAwAkECABiBIAEAjECQAABGIEgAACMQJACAEQgSAMAI\nBAkAYASCBAAwAkECABiBIAEAjECQAABGIEgAACMQJACAEQgSAMAI/x+apbG7WG4p5QAAAABJRU5E\nrkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "% We switch the local indices of nodal star to the first case. \n",
    "ix = (elem(nodeStar(1,:),3) == elem(nodeStar(4,:),2));\n",
    "iy = (elem(nodeStar(2,:),2) == elem(nodeStar(3,:),3));\n",
    "% nodeStar(:,ix & iy)  = nodeStar([1 2 3 4],ix & iy);\n",
    "nodeStar(:,ix & ~iy) = nodeStar([1 3 2 4],ix & ~iy);\n",
    "nodeStar(:,~ix & iy) = nodeStar([1 4 2 3],~ix & iy);\n",
    "% nodeStar(:,~ix & ~iy)= nodeStar([1 4 3 2],~ix & ~iy);\n",
    "%%\n",
    "% We label the local indices of nodes in the star in the following figure.\n",
    "figure(2)\n",
    "node1 = [0,0; -1,0; 1,0; 0,1; 0,-1];\n",
    "elem1 = [1 4 2; 1 3 4; 1 2 5; 1 5 3];\n",
    "showmesh(node1,elem1); findelem(node1,elem1); findnode(node1);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then t1 and t2 are merged to a big triangle and stored in t1 and\n",
    "similarly t3 and t4 are merged together and stored in t3. Elements t2 and t4 will be removed later on."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "t1 = nodeStar(1,:); \n",
    "t2 = nodeStar(2,:); \n",
    "t3 = nodeStar(3,:);\n",
    "t4 = nodeStar(4,:);\n",
    "% p1 = elem(t1,1); \n",
    "p2 = elem(t1,3); \n",
    "p3 = elem(t2,2); \n",
    "p4 = elem(t1,2); \n",
    "p5 = elem(t3,2);\n",
    "elem(t1,:) = [p4 p2 p3]; \n",
    "elem(t2,1) = 0;\n",
    "elem(t3,:) = [p5 p3 p2]; \n",
    "elem(t4,1) = 0;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For boundary good nodes, there is only one case."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA\nB3RJTUUH4QgJBy8n7u39TwAAACR0RVh0U29mdHdhcmUATUFUTEFCLCBUaGUgTWF0aFdvcmtzLCBJ\nbmMuPFjdGAAAACJ0RVh0Q3JlYXRpb24gVGltZQAwOS1BdWctMjAxNyAxNTo0NzozOROhAJoAABcx\nSURBVHic7d1/jJz1feDx78YbbA5jHKqgHPaGNViNFFeqObUBnIDXJL4ktIAJjaE0EG+IExRVkYnp\nH1WiYDenqn/EGHKhoigCkx89xVwApxINIIVxICxGcc7Q+iAJvky09R3BOKzXLLgUe+6PxyzL/vbu\nzDyf53leL/mP8ezj4cN+Z5/3PjPPzHQ0Go0EAHl7R94DAEBKggRAEIIEQAiCBEAIggRACIIEQAiC\nBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIE\nQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRA\nCIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAI\nggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiC\nBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIEQAiCBEAIggRACIIE\nQAiCBEAIggRACIIEQAiCBEAIggRACIIEebr77rsvvfTSSy+99JFHHsl7FshZZ94DQHXt3r37hhtu\neP3111NKa9asyXscyJkjJMjH4ODgVVddldUISIIEeVm/fv2+ffvyngICESTIwR133LF9+/Z58+Zd\nd911ec8CUXgOCdrt6aefvvHGG1NKt9122yuvvJL3OBCFIyRorXq9PvKvr7zyytq1a48cOXL11Vd/\n7nOfy2koiEiQoFXq9fqqVauWLFnS0dGxatWqrEw33HDDL3/5y6VLl9555515DwixeMgOWqW3t7dW\nq2WXa7Vab2/vtdde+73vfW/u3Lnbt28/9dRTc50OwuloNBp5zwDl1NHRMeqaT37yk/fee+/JJ5+8\naNGi7JqXX3754MGDKaUzzjhjwYIFt9122yWXXNLuQSEGR0jQEqOeOsoMDQ2llF577bXnn39+1Jde\nfPHFF1988fDhw22YDWISJGiJ7u7unp6e4YfsUko9PT1r1qxZvHjxyM2eeeaZJ598MqW0cuXK973v\nfeecc06b54Q4BAlaol6v79n38+Vrl+3Zvjel1H1B18DR361evXr9+vUjN7vllluyIF177bXXX399\nPrNCDIIEzVev18+96A/XbP1Y94quhV0LUko9G1cM9A9ecd3l9397R3d3d94DQkRO+4YmG1mjkdcv\n7FrQ8/Xzrrju8nGfXgKcZQfNNLZGtS1PpJR6Nq7I/jrQP1i7aZfjJBjLERI0zUTHRiM5ToKJCBI0\nx3RqlNEkGJcgQRNMv0YZTYKxBAlm60RrlNEkGEWQYFZmVqOMJsFIggQzN5saZTQJhgkSzNDsa5TR\nJMgIEsxEs2qU0SRIggQz0NwaZTQJBAlOTCtqlNEkKk6Q4AS0rkYZTaLKBAmmq9U1ymgSlSVIMC3t\nqVFGk6gmQYKptbNGGU2iggQJptD+GmU0iaoRJJhMXjXKaBKVIkgwoXxrlNEkqkOQYHwRapTRJCpC\nkGAccWqU0SSqQJBgtGg1ymgSpSdI8DYxa5TRJMpNkOAtkWuU0SRKTJDguPg1ymgSZSVIkFJxapTR\nJEpJkKBgNcpoEuUjSFRdEWuU0SRKRpCotOLWKKNJlIkgUV1Fr1FGkygNQaKiylGjjCZRDoJEFZWp\nRhlNogQEicopX40ymkTRCRLVUtYaZTSJQhMkKqTcNcpoEsUlSFRFFWqU0SQKSpCohOrUKKNJFJEg\nUX5Vq1FGkygcQaLkqlmjjCZRLIJEmVW5RhlNokAEidJSo4wmURSCRDmp0UiaRCEIEiWkRmNpEvEJ\nEmWjRhPRJIITJEpFjSanSUQmSJSHGk2HJhGWIFESajR9mkRMgkQZqNGJ0iQCEiQKT41mRpOIRpAo\nNjWaDU0iFEGiwNRo9jSJOASJolKjZtEkghAkCkmNmkuTiECQKB41agVNIneCRMGoUetoEvkSJIpE\njVpNk8iRIFEYatQemkReBIliUKN20iRyIUgUgBq1nybRfoJEdGqUF02izQSJ0NQoX5pEOwkScalR\nBJpE2wgSQalRHJpEewgSEalRNJpEGwgS4ahRTJpEqwkSsahRZJpESwkSgahRfJpE6wgSUahRUWgS\nLSJIhKBGxaJJtIIgkT81KiJNoukEiZypUXFpEs0lSORJjYpOk2giQSI3alQOmkSzCBL5UKMy0SSa\nQpDIgRqVjyYxe4JEu6lRWWkSsyRItJUalZsmMRuCRPuoURVoEjMmSLSJGlWHJjEzgkQ7qFHVaBIz\nIEi0nBpVkyZxogSJ1lKjKtMkTogg0UJqhCYxfYJEq6gRGU1imgSJllAjRtIkpkOQaD41YixNYkqC\nRJOpERPRJCYnSDSTGjE5TWISgkTTqBHToUlMRJBoDjVi+jSJcQkSTaBGnChNYixBYrbUiJnRJEYR\nJGZFjZgNTWIkQWLm1IjZ0ySGCRIzpEY0iyaRESRmQo1oLk0iCRIzoEa0giYhSJwYNaJ1NKniBIkT\noEa0miZVmSAxXWpEe2hSZQkS06JGtJMmVZMgMTU1ov00qYIEiSmoEXnRpKoRJCajRuRLkypFkJiQ\nGhGBJlWHIDE+NSIOTaoIQWIcakQ0mlQFgsRoakRMmlR6gsTbqBGRaVK5CRJvUSPi06QSEySOUyOK\nQpPKSpBISY0oGk0qJUFCjSgkTSofQao6NaK4NKlkBKnS1Iii06QyEaTqUiPKQZNKQ5AqSo0oE00q\nB0GqIjWifDSpBASpctSIstKkohOkalEjyk2TCk2QKkSNqAJNKi5Bqgo1ojo0qaAEqRLUiKrRpCIS\npPJTI6pJkwpHkEpOjagyTSoWQSozNQJNKpDOvAegVYpSo4H+wT3b/3X4cveKruVrl+U7EpMr3JIN\nN+n+b+/o7u7Oexwm1NFoNPKegeYLXqOB/sGB/kP1vv492/cO9A+O3WBh14LuC7qWr10Wc/4TUtvy\nREqpZ+OKvAeZlRIs2UD/YO2mXZoUmSCVUPwa1bY8sWf73uFrurvTunVp5cq0c2eq1VKt9tbG3Rd0\nrbn1Ywu7FrR/zmYpQZBKs2SaFJwglU38Gj2w4Uf1vv7hPVpPz+ht6vVUq6WdO49fWNi1YN3/vCrm\nDm46ih6kki2ZJkUmSKVSlBqtW5fuvnta/2Tz5rRpU+gd3JQKHaRSLpkmheUsu/IoX41SSjffnDZt\nSgP9g9v+7PvjPnVB65R1yZx3F5YglUQpa5SJv4MrpXIvmSbFJEhlELxGKaXalidmtmvLDO/gsoe/\naIPSL5kmBSRIhRe/Rnu2781O0JrZri1z882ppyft2b63/kR/0yZjAhVZMk2KRpCKLX6NUkrZ/mjd\nuik2++u/TitWpKeemnCD7OSuel/QvVuZTLlk99yTVqxICxakpUvTX/1VOnRo/M3iL5kmhSJIBVaI\nGqU390crV062zVNPpdtuS319E+7ahm9h5KthaJHJl+xv/zatW5f6+tLhw2nfvvT1r6fLLktHj46z\nZSGWTJPiEKSiKkqNUkrZ09oT/bp9663pyivThRem116b4nayc3QH+gdjPk9eJpMs2e9+l/7u71JK\nadOmNDCQvvGNlFL6yU/Sj388zsZFWTJNCkKQCqlANcp+O57kwZ9t29J996XXX5/6prq7jz8ENNA/\n8WEUszb5kv385+nw4XT66enLX06nnZb+8i/T4sUppfQv/zLOxgVaMk2KQJCKp0A1Sm8+GzHJCxBv\nuSXdd1+6775p3Vq2d6tt6WvCZExg8iXr6Egf/3j6i79InZ0ppfTGG2loKKWUFi0af/sCLZkm5c67\nfRdMsWo07KyzJvzSxRef8O3EfP1/yUy0ZB/+cPrwh49fbjTShg3p5ZfTwoUTrmOxlsz7gufLEVKR\nFLFG2ZMHv/lNc26tWbfDJKa5ZAcOpD/90/T3f59OOinddVd697vH36xwS+Y4KUeCVBhFrFFKqWfj\nBSm97d2gZyO7nWJ9BwpnOkvW15fOPTc9+GBavDg98ki64ooJtyzikmlSXgSpGApao5TSwq7TUkrN\n+rk+vne7oGDfhGKZcskefTR95CNp//708Y+np59OF1002a0VdMk0KReCVADFrVF683Pbsk8lmKVt\n247fYFGekCioyZdsaChdc0169dX0J3+S/umf0umnT3ZThV4yTWo/QYqu0DXKdK9YnFK6557Z3s7O\nnSkV8HftIppkyb773fTCCyml9M//nObOTZ2dx/9s3DjOxkVfMk1qM0EKrQQ1Sm/uj2Z/hFTEZyMK\napIle+aZ4xeOHUtHj77159ixcTYuwZJpUjsJUlzlqFFKaWHXacvXLqvX0+bNk23WaKRGI61ePf5X\ne3tTvZ6Wr122fO2yVgzJSJMs2e23H1+pUX+2bh29ZWmWTJPaRpCCKk2NUvbzvHHFwq4FmzZN0aSJ\n9PambdvSwq4Fa279WLOnYxyWbBRNag9BiqhMNcoMf6D1DHZww7u2DbvWt2Y6xmHJRtGkNhCkcMpX\no8zMdnCl3LUVhSUbRZNaraPRaOQ9A28pa42GDX+sdfa2mytXjvMmnvV6qtfTzp1p27ZUrxd+15Z9\nZGrPxhV5DzJDFVyyyQ30D9Zu2uW9hVpBkAIpfY0yA/2DD2z40fCHtmU/1OvWpZUr086dqVZ728ld\n3Rd0rfvB2vYP2URFD1Kq3pJNSZNaRJCiqEiNMgP9gwP9h/Zs3zvQPzjux4ku7FqwfO2y7gu6SvDd\nKEGQUsWWbDo0qRUEKYRK1WiUbAeXfeRB9nr+ou+7RylHkEYq/ZJNkyY1nSDlr8o1qoLyBYlhmtRc\nzrLLmRpBcTnvrrkEKU9qBEWnSU0kSLlRIygHTWoWQcqHGkGZaFJTCFIO1AjKR5NmT5DaTY2grDRp\nlgSprdQIyk2TZkOQ2keNoAo0acYEqU3UCKpDk2ZGkNpBjaBqNGkGBKnl1AiqSZNOlCC1lhpBlWnS\nCRGkFlIjQJOmT5BaRY2AjCZNkyC1hBoBI2nSdAhS86kRMJYmTUmQmkyNgIlo0uQEqZnUCJicJk1C\nkJpGjYDp0KSJCFJzqBEwfZo0LkFqAjUCTpQmjSVIs6VGwMxo0iiCNCtqBMyGJo0kSDOnRsDsadIw\nQZohNQKaRZMygjQTagQ0lyYlQZoBNQJaQZME6cSoEdA6FW+SIJ0ANQJarcpNEqTpUiOgPSrbJEGa\nFjUC2qmaTRKkqakR0H4VbJIgTUGNgLxUrUmCNBk1AvJVqSYJ0oTUCIigOk0SpPGpERBHRZokSONQ\nIyCaKjRJkEZTIyCm0jdJkN5GjYDIyt0kQXqLGgHxlbhJgnScGgFFUdYmCVJKagQUTSmbJEhqBBRS\n+ZpU9SCpEVBcJWtSpYOkRkDRlalJ1Q2SGgHlUJomVTRIagSUSTmaVK0g/eAHP/jEJz6xePHis885\n+6RTTjrwy4PH3jiW91BU1wt7X/zHT9//j5++/8Cvfpf3LBTeDJp09OjRrVu3XnTRRaeffvrixYsv\nueSShx9+uJUzTqGj0Wjk+J9vm0aj8YUvfOGOO+4Ydf2i5e/pvf+qzrmduUxFFdS2PJFS6tm4YuyX\ndtz40P/6/r+mlD6z4+r3/vGidk9GGQ30D9Zu2nX/t3d0d3dPvuWRI0cuvvjivr6+Uddv2LBh69at\nrZpvUlU5Qrr33nuP16gjdf3RmV1/dGZ2/f49Lzz8NzvznIxKGjr42sN/szOrETTR9I+Tvva1r2U1\nmjdvXm9v72WXXXbSSSellG699dYHHnigDaOOVZUg3XzzzdmFVTetuP6Hf379D//8E9+8JLtm7w9/\nkd9cVM4z9z373z941y3/5Y4n7vhZ3rNQTtNpUqPR+OY3v5lS6uzs3Ldv31133bVjx46vfvWr2Vcf\nfPDB9ow6SiWCNDQ09Nxzz6WUFvznU1feeEF25bJLf3/egrkppaGDrx3+7VCe81ElQwdePfjrl4/+\nhycvaaEpm/SrX/1qcHAwpXThhReeeebxR4wuv/zy7MKRI0faMuZolXjuZNeuXR0dHZ3z5px1/uLh\nKxvHGm/8+9GUUufczlN+7+T8pqNa3nveoo98+cLs8u7vPvPybw7lOw9lNdykcZ9PWrJkyW9/+9uU\n0sknH9/7HTt27M4778wur169uo2TvqW0QarX69ka1Ov1K9dd8el7PznqDO/HvrHrjX9/I6XU9cdn\nvqOzEkeKRLBo+XsWLX9Pdvn5H/9akGidUU0a3iumlN75zneeccYZ2eVnn33285///HPPPXfgwIE5\nc+Z85Stfueaaa3IZuIQ74nq9vmrVqiVLlnR0dKxatWrs642O/sfRhzbVdm59MqXUedKcS/7bxfkN\nS5kN9A9uu3J7bUtfbUvftiu3D/QP5j0RlTPcpJF7xVGP4x0+fPixxx47cOBASmn+/PlnnXXWO96R\nTxpKGKTe3t5arZZdrtVqCxefNrJGL/zvA//w0e/23bk7pXTywnmf+h9/9u7f/71c5qT0Htjwo3pf\nf3a53tf/wIYf5TsP1bSwa8HAGy+P3Cv29vaO3GDp0qXf+c53Nm/e3NXVdejQoc985jPr16/PYdBS\nvg6po6Nj1DU9G4+fyPBvP/9//+cnvzl2tJFSetdZp73vv56TndcArVDbMvoVHsN3xcye7+8d+LfB\nlNK5V//BaYtObd9kVMzYu+Krr776+uuvp5Tmz58/Z86c7Mr9+/d3dXU1Go1TTjnlpZdemjdvXpvn\nLNtzSOOeUtJz6kcbjcZDDz30/JP1lFJnZ+fq1avPO++8Ns9G1dTSmCCd+tGRf/31nP87kAZTSuf+\npw+899T3tm8yKmbsXfFLX/pS9tLMHTt2XHbZZdmVixYtOvPMM/fv3z80NLR79+4PfvCD7R60UTo9\nPT0j/wd7enoajcbwC49POeWUn/70p3nPSCWMe1ccaeXKldmXHn/88VwmpCLG3hXvueee7PJ11103\nvNlTTz01vM3BgwfbP2cJH7Kr1+ubN2/etm1bSqmnp+fRRx996aWXli5deujQoZTS+eeff/7554/6\nJ1/84heXLFnS/lEpt7F3xVEb9PT07Ny5M6X0+OOP5/DbKJUx9q64f//+97///dlLkTZs2HDllVc+\n/fTTt99++7PPPptS+sAHPrBr1672z1nCIGVGnuC4Y8eONWvWTLLxY4899qEPfagdY1E9I++KowgS\n7TTqrvitb31r3JMX3vWud/3sZz87++yz2zfZm0p4ll1m5Pd97969+Q1C1U35HpfQHqPuip/97Gcf\neeSRZcuWDV8zZ86cT33qU7t3786lRqnER0gATMfAwMAvfvGL+fPnL126dO7cPE88FiQAQijtQ3YA\nFIsgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQ\ngiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCC\nIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIg\nARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiAB\nEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQ\ngiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCC\nIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIgARCCIAEQgiABEIIg\nARCCIAEQwv8HFbOUHDF/IvUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "figure(2)\n",
    "node1 = [0,0; -1,0; 1,0; 0,1];\n",
    "elem1 = [1 4 2; 1 3 4];\n",
    "showmesh(node1,elem1); findelem(node1,elem1); findnode(node1);\n",
    "[ii,jj] = find(t2v(:,bdGoodNode));\n",
    "nodeStar = reshape(ii,2,length(bdGoodNode));\n",
    "t1 = nodeStar(1,:); \n",
    "t2 = nodeStar(2,:);\n",
    "p1 = elem(t1,1); \n",
    "p2 = elem(t1,3); \n",
    "p3 = elem(t2,2); \n",
    "p4 = elem(t1,2);\n",
    "elem(t1,:) = [p4 p2 p3]; \n",
    "elem(t2,1) = 0;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clean and Shift Node and Elem Matrices\n",
    "\n",
    "The empty rows in node and elem matrices should be relased for the\n",
    "efficent usage of meomory. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "elem((elem(:,1) == 0),:) = [];\n",
    "isGoodNode = false(N,1);\n",
    "isGoodNode(intGoodNode) = true;\n",
    "isGoodNode(bdGoodNode) = true;\n",
    "node(isGoodNode,:) = [];"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the indices of nodes in the coarsened mesh is changed while\n",
    "elem still use the indices in the fine mesh. We build the indices map\n",
    "between the fine grid to the coarse grid. For example, `indexMap(10) = 6`\n",
    "means the 10-th node in the fine grid  is now the 6-th node in the coarse\n",
    "one."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "indexMap = zeros(N,1);\n",
    "indexMap(~isGoodNode)= 1:size(node,1);\n",
    "elem = indexMap(elem);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the application to adaptive finite element methods, we will remove\n",
    " good-to-coraen nodes whose star are marked for coarsening. \n",
    "\n",
    "In the application to multigrid methods, we will record all neighboring\n",
    "nodes in the nodal star of good nodes; see the functoin: `uniformcoarsen`.\n",
    "\n",
    "The 3-D coarsening algorithm is slightly complicated than the 2-D case."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aditional data structure\n",
    "\n",
    "- `bdFlag` stores information on boundary conditions. It is updated as"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```bdFlag(t1,:) = [bdFlag(t1,2) bdFlag(t2,1) bdFlag(t1,1)];\n",
    "bdFlag(t3,:) = [bdFlag(t3,2) bdFlag(t4,1) bdFlag(t3,1)];\n",
    "bdFlag(t5,:) = [bdFlag(t5,2) bdFlag(t6,1) bdFlag(t5,1)];\n",
    "bdFlag((elem(:,1) == 0),:) = [];```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- `tree(:,1:3)` stores the binary tree of the coarsening. `tree(:,1)` is the index of parent element in the coarsened mesh and `tree(:,2:3)` are two children indices in the original mesh."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```tree = zeros(NTdead,3,'uint32');\n",
    "tree(1:NTdead,1) = [t1'; t3'; t5(idx)'];\n",
    "tree(1:NTdead,2) = [t1'; t3'; t5(idx)'];\n",
    "tree(1:NTdead,3) = [t2'; t4'; t6(idx)'];```"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Matlab",
   "language": "matlab",
   "name": "matlab"
  },
  "language_info": {
   "codemirror_mode": "octave",
   "file_extension": ".m",
   "help_links": [
    {
     "text": "MetaKernel Magics",
     "url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md"
    }
   ],
   "mimetype": "text/x-matlab",
   "name": "matlab",
   "version": "0.14.3"
  },
  "toc": {
   "colors": {
    "hover_highlight": "#DAA520",
    "navigate_num": "#000000",
    "navigate_text": "#333333",
    "running_highlight": "#FF0000",
    "selected_highlight": "#FFD700",
    "sidebar_border": "#EEEEEE",
    "wrapper_background": "#FFFFFF"
   },
   "moveMenuLeft": true,
   "nav_menu": {
    "height": "138px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": 4,
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": true,
   "widenNotebook": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
